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Laminar  flow  of  fluid  through  fibrous  and  granular  media  and  deposition  of 
colloidal  particles  from  a liquid  suspension  are  two  fundamental  phenomena  encountered 
in  many  industrial  applications.  Finite  Element  Method  (FEM)  and  Effective  Medium 
Approximation  (EMA)  are  used  to  determine  the  fluid  flow  permeability  and  particle 
capture  efficiency  of  regular  array  and  random  arrays  of  cylindrical  and  spherical 
collectors.  The  FEM  takes  into  account  of  the  detailed  layout  of  the  packing  but  it  can 
not  handle  randomly  packed  collectors  yet  at  the  current  stage,  while  the  EMA  gives  an 
averaged  estimation  of  all  the  flow  parameters  but  is  not  able  to  provide  a resolution  to 
distinguish  the  detailed  structure  differences  and  orientations. 

The  EMA  assumes  a model  system  in  which  a packing  element  (a  single  fiber  in 
the  fibrous  medium  and  a single  sphere  in  the  granular  medium)  is  surrounded  by  a fluid 


xv 


envelope  and  an  effective-medium  beyond  the  envelope.  It  integrates  the  important 
features  of  both  the  cell  models  and  Brinkman’s  model.  The  Stokes  equation  and 
Brinkman’s  equation  are  solved  for  the  fluid  envelope  and  effective  medium  regions, 
respectively,  to  obtain  the  permeability  and  close-to-surface  velocity  field  around  the 
collectors.  The  convective  diffusion  equation  is  then  solved  to  determine  the  particle 
deposition  rate.  The  analytical  expressions  for  the  permeability  and  particle  deposition 
rate  are  derived  for  all  possible  cases  of  random  packing  of  uniform  and  non-uniform 
cylinders  and  spheres. 

Effects  of  various  system  properties  and  operating  conditions  on  deposition  of 
colloidal  particles  are  investigated.  The  physical  or  chemical  conditions  include  the 
properties  which  affect  the  magnitude  of  double  layer  interaction,  the  electrolyte 
concentration  and  surface  potentials,  and  the  property  which  affects  the  van  der  Waals 
interaction,  the  Hamaker  constant.  It  was  found  that  the  effects  of  the  above  properties 
are  much  more  significant  when  the  surface  interactions  play  more  important  roles  in  the 
particle  deposition,  i.e.,  when  the  height  of  the  total  interaction  energy  barrier  is  higher 
than  5 kBT.  Particle  deposition  becomes  virtually  impossible  when  the  height  of  the 
repulsive  energy  barrier  increases  beyond  20  k0T. 
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CHAPTER  1 

GENERAL  INTRODUCTION 


1-1  Fibrous  and  Granular  Porous  Medium  Models 


Theory  of  fluid  flow  through  fibrous  and  granular  media  is  a branch  of  the  physics 
of  fluid  flow  through  porous  media.  In  the  current  study,  our  major  interest  is  the  slow 
flow  with  applications  to  particle  deposition  in  filter  media,  although  the  developed 
results  can  have  much  broader  uses.  In  this  instance,  the  concept  of  laminar, 
incompressible  and  Newtonian  flow  though  porous  media  governed  macroscopically  by 
Darcy’s  law  can  be  applied.  The  Darcy’s  law  simply  states  that  the  pressure  gradient 
across  a porous  medium  is  proportional  to  the  superficial  velocity  and  the  fluid  viscosity: 


— = — Uu 

dl  k 


(1-1) 


the  reciprocal  of  the  proportionality  constant,  k,  is  called  the  permeability  of  the  medium. 
On  the  other  hand,  the  Stokes  law  applies  in  the  clear  fluid  regions  locally: 

Vp  = pV2u  (1-2) 

The  solution  of  the  Stokes  equation  (1-2)  gives  the  fluid  velocity  field  needed  for  many 
mass  and  heat  transfer  applications.  The  difficulty  in  applying  equation  (1-2)  in  the 
fibrous  and  granular  media  is  with  the  boundary  conditions.  Because  of  the  complex 
geometric  nature  of  the  fluid-solid  boundaries  in  the  media,  some  simplified  models  have 
to  be  introduced.  There  have  been  a few  good  reviews  on  the  porous  medium  models. 


1 
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such  as  Rajagopalan  and  Tien  (1979),  Tien  (1989)  and  Pich  (1986).  Summarized  below 
are  the  key  concepts  and  the  major  conclusions  of  the  models  which  have  applications  in 
mass  and  heat  transfer  processes  and  which  will  be  referred  to  frequently  in  the  context  of 
the  current  work  for  comparison. 

The  early  porous  medium  models,  such  as  the  well  known  Kozeny-Carman’s 
model  based  on  the  equivalent  capillary  idea,  concerned  only  the  pressure  gradient 
predictions.  The  models  developed  later  on,  like  the  Happel’s  and  Kuwabara’s,  can  not 
only  give  reasonable  good  predictions  of  pressure  gradient  but  also  a representative  flow 
field  in  the  vicinity  of  the  individual  elements.  It  is  this  second  class  of  models  that  is  of 
interest  to  the  current  study. 

The  widely  used  cell  models  portray  a porous  medium  as  a collection  of  unit  cells 
of  a finite  size  which  consist  of  a spherical  or  cylindrical  solid  element  surrounded  by  a 
fluid  envelope.  The  size  of  the  solid  element  is  assumed  to  be  identical  to  the 
characteristic  size  of  the  packing  elements,  and  the  size  of  the  fluid  envelope  is 
determined  based  on  the  porosity  of  the  porous  medium  (Happel  1958).  For  the 
description  of  the  flow  characteristics  in  a porous  medium,  the  momentum  equation  is 
solved  for  the  simplified  flow  geometry  with  appropriate  boundary  conditions.  Since  the 
characteristic  pore  size  of  packed  beds  are  usually  small,  the  Reynolds  number  based  on 
the  pore  size  is  small.  Consequently,  the  creeping  flow  equation  is  typically  solved  in 
these  models  with  the  no-slip  condition  at  the  solid  surface  (i.e.,  at  the  inside  wall  of 
capillary  models  or  at  the  surface  of  the  solid  element  of  cell  models).  Another  boundary 
condition  needs  to  be  specified  at  the  outer  surface  of  the  fluid  envelope.  In  Happel’s 
model  (1958),  a stress-free  condition  is  applied  whereas  a vorticity-free  condition  is 


3 


applied  in  Kuwabara’s  model  (1959).  Although  either  conditions  may  appear  to  be 
plausible  considering  the  symmetry  at  the  surface  of  the  fluid  envelope,  some  ambiguity 
exists  regarding  the  boundary  conditions  and  these  cell  models  may  not  account  for  the 
influence  of  the  neighboring  elements. 

The  physical  picture  underlying  the  Brinkman’s  model  (Brinkman  1947)  is  that  of 
a packing  element  embedded  in  a macroscopic  porous  medium.  The  flow  field  is 
described  by  the  Darcy’s  law  at  distances  remote  from  the  element,  but  by  the  Stokes 
equation  near  the  element;  in  other  words,  the  Brinkman’s  equation: 


is  solved.  The  essence  of  the  Brinkman  equation  is  that,  on  the  average,  the  fluid  in  the 
proximity  to  an  obstacle  embedded  in  a porous  medium  experiences  a body  damping 


accounts  for  the  influence  of  neighboring  solid  elements  on  the  flow.  The  Brinkman’s 
model  applied  to  fibrous  media  was  worked  out  by  Spielman  and  Goren  (1968),  to 
Granular  media  by  Payatakes  et  al.  (1974b). 

The  velocity  fields  in  the  vicinity  of  the  fiber  surfaces  for  all  these  models  in  the 
case  of  flow  perpendicular  to  all  the  fiber  axes  in  fibrous  media  can  be  summarized  by  the 
following  equation: 


1 2 

Vp  = pu  + pV  u 


k 


(1-3) 


force  proportional  to  the  velocity  in  addition  to  the  viscous  and  the  pressure  forces,  which 


H/(r,0)  = 


aU  sin  0 a 


(1-4) 


2A  r 


r a a a 
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where  r and  0 are  the  polar  coordinates  and  vj/(r,0)  is  the  stream  function  related  to 
velocity  components  by  u r = (1  / r)(3i|/  / 30)  and  ue  = -3t|/  / 3r , a is  the  radius  of  the 
packing  elements,  A is  the  model  parameter  given  by 

A = -^-^ln(l  — e)  + (1  — e)  — ^(1  — e)2  (Kuwabara’s  model)  (1-5) 


a 1 , (1-s) 

A = ln(l-e)  + 


A = • 


2 2 l + (l-e)' 

K0(a/Vk) 


(Happel’s  model) 


(Brinkman’s  model) 


(1-6) 


(1-7) 


(a/Vk)K,(a/Vk) 

For  some  applications  it  is  adequate  to  approximate  the  velocity  field  of  equation  (1-4)  by 
using  only  the  first  term  in  a power  series  expansion  by  powers  (r-a)/a  . The  result  is 

-|2 

io  r 
A a 


M/(r»0)  = 


aU  sin0 


(1-8) 


When  the  flow  is  parallel  to  the  fibers,  Happel’s  model  gives  the  following  flow  field 


u(r)  = - 


2U 

2 2 2a2  r 

a - r + In- 

1-s  a 

a2(l-e) 

2 In  1 3 + 4(1  e) 

1 - 8 

(1-8)2 

(1-9) 


The  flow  field  predicted  by  the  Brinkman’s  model  for  the  current  flow  pattern  has  not 
been  reported,  but  the  result  will  be  given  in  Chapter  3 considering  the  fact  that  it  is  a 
special  case  of  the  effective  medium  approximation  model  studied  in  the  current  research. 
The  permeability  predictions  by  these  models  are  summarized  in  Table  3-1  and 


Table  3-2. 
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The  general  velocity  field  around  the  spheres  in  granular  media  for  the  cell 
models  can  be  written  as: 


„ a r 

(A2  (A 

4' 

v = 

Ci-  + C2-+C3 

- +c4  - 

r a 

W \aJ 

sin2  0 


(1-10) 


where  the  stream  function  \\i  is  defined  by  ur  = [-1  / (r  sin0)](3vy  / 30)  and 
ue  = [1  / (rsin0)](3\|r  / 3r) , and  Ch  C2,  C3  and  C4  are  the  model  parameters.  For 
Happel’s  model, 


C,  = 


Ua 


2 


2 - 3(1  - e)1/3  + 3(1  - s)5/3  - 2(1  - s)2 


(1-11) 


Co  = 


Ua" 


3 + 2(l-s) 


5/3 


2 - 3(1  - s)1/3  + 3(1  - s)J/J  - 2(1  - e) 


5/3 


-1 


(1-12) 


C,= 


Ua" 


2 + 3(1  - s)5/3  2 - 3(1  - e)1,J  + 3(1  - s)3,J  - 2(1  - e) 


1/3 


v5/3 


-1 


(1-13) 


C4  =-^—(1  -e)5/3 


and  for  Kuwabara’s  model. 


C,  = 


2 - 3(1  - s)1/3  + 3(1  - s)5/3  - 2(1  - e)2 


Ua2[3  + 2s] 


20 


2 - 3(1  - s)1/3  + 3(1  - e)5/3  - 2(1  - e)2 


l-d-E) 


1/3 


(1-14) 


(1-15) 


C2  = 


- 3Ua" 


2 - 3(1  - e)i/3  + 3(1  - e)5/3  —2(1  — E)2]fl  — (1  — e)1/3 


Ci  = 


Ua2[l  + e] 


2 - 3(1  - e)1/3  + 3(1  - s)5/3  - 2(1  - e)2 


l-(l-s) 


1/3 


(1-16) 


(1-17) 
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CA  = 


20 


- 3Ua^(l  - s) 

2 - 3(1  - s)1/3  + 3(1  - s)5/3  - 2(1  - s)2  lfl  — (1  — e)1/3 


i3 


(1-18) 


The  stream  function  for  the  Brinkman’s  model  in  granular  media  is  (Payatakes  et  al. 
1974b): 


DO 

1 

II 

> 

0 + 2m3 

fl+m3m?  e_m,r(l  + m,r)l 

z. 

m,  r 

l m2  ) 

sin2  0 


(1-19) 


where 


m,  = 


Vk 


(1-20) 


m2  — em' 


3 

2 


(1-21) 


m3  = — 
3 2 


-l  + 4(emi-l-m1) 

mt  ' ’ 


(1-22) 


The  permeability  predictions  of  flow  through  granular  media  are: 
.2 


2a2 


k = 


1 - 3(1  - s)1/3  + 3(1  - e)5/3  + 2(1  - s)2 


3 3 + 2(l-s) 


5/3 


(Happel’s  Model)  (1-23) 


k = 


2 a 
9 1-s 


B(8) 


(Kuwabara’s  Model)  (1-24) 


where 


B(8)  = 


l-(l-8) 


1/3 


l + |(l-8)'/3+|(l-8)2/3+i(l-8) 


-i2 


1 + 11(1  ■ - e)1'3  + T(1 . e)M  _ JL(1 . £) +2(1  _ s)« 

50  50  25  5 

+ li(1_e)2+A(l_s)’«+i_(1_s)8/3 
L 25  25  25 


(1-25) 
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Theoretically  speaking,  the  cell  models  and  the  Brinkman’s  model  can  be 
employed  only  in  cases  where  the  packing  porosity  is  relatively  high  (e.g.,  s>0.6,  Happel 
and  Brenner  1965,  Spielman  and  Goren  1968,  Pich  1986)  although  all  of  them  have  been 
used  in  the  study  of  granular  beds  where  the  porosity  is  usually  in  the  range  of  0.35-0.45, 
simply  because  no  other  good  model  is  available  for  low  porosity  packings  where  the 
effect  of  the  neighboring  elements  plays  a very  important  role.  The  theoretical  models 
for  packings  of  non-uniform  fibers  and  spheres  have  not  been  found  in  the  literature. 


1-2  Particle-Collector  Interactions 

Dispersion  of  particles  in  liquids  is  commonly  encountered  in  a wide  range  of 
process  industries.  Dispersed  particles  can  have  sizes  ranging  from  fractions  of  a 
millimeter  to  macromolecular  dimensions  (a  few  nanometers).  Particles  bigger  than 
about  1pm  in  size  are  usually  called  suspended  particles,  smaller  than  about  lnm  are 
molecules  or  atoms  while  the  sizes  in  between  are  usually  defined  as  colloidal  particles. 

Compared  with  suspended  particles,  colloidal  particles  may  have  the  following 
properties.  They  are  difficult  to  be  seen  under  the  usual  optical  microscopes  because  their 
sizes  are  mostly  smaller  than  the  wavelength  of  the  ordinary  light  (about  0.5pm). 
Sedimentation  can  not  be  used  to  separate  them  from  the  suspending  liquids  because  of 
the  unacceptable  long  time  required.  It  can  be  shown  that  the  important  types  of  inter- 
particle (or  particle-surface)  forces,  known  collectively  as  colloidal  interactions,  vary 
roughly  in  proportion  to  the  particle  size,  whereas  external  forces,  such  as  gravity  and 
hydrodynamic  drag  have  a stronger  dependence  on  particle  size.  For  instance,  the  gravity 
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force  is  proportional  to  the  mass  of  particle  and  hence  to  the  cube  of  the  particle  size.  The 
hydrodynamic  forces,  as  a result  of  fluid  flow,  depend  roughly  on  the  square  of  particle 
size.  This  means  that  for  large  particles,  gravitational  and  drag  forces  will  predominate 
over  colloidal  forces,  but  the  later  may  be  much  more  significant  for  colloidal  particles. 
Colloidal  interactions  play  a crucial  role  in  determining  the  particle  deposition  efficiency, 
stability  and  rheological  behavior  of  colloidal  dispersions. 

There  are  several  types  of  colloidal  interactions  which  are  known  to  be  important, 
the  best  known  of  which  are  the  van  der  Waals  attraction  and  double  layer  interactions. 
These  two  contributions  form  the  basis  of  the  classical  DLVO  theory  of  colloidal 
stability,  independently  proposed  by  Derjaguin  and  Landau  (1941)  and  Verwey  and 
Overbeek  (1948).  Other  effects  are  now  known  to  be  important  in  certain  cases  and  may 
play  a large  part  in  determining  whether  or  not  particles  will  adhere  to  each  other  or  to 
other  surfaces.  In  aqueous  systems,  the  effect  of  hydration  can  be  important.  This  is  often 
associated  with  hydration  of  ions  at  particle  surfaces  and  usually  gives  an  extra  repulsion. 
It  is  now  known  that  hydrophobic  effects  can  also  be  important,  giving  an  extra  attraction 
between  particles  or  between  particle  and  other  surfaces.  Other  significant  effects  arise 
from  the  presence  of  adsorbed  polymers,  giving  either  a repulsion  (steric  interaction)  or 
an  attraction  (polymer  bridging).  There  have  been  a few  very  good  reviews  on  this 
subjects  (Isrealachvili  1992,  Russel  et  al.  1989).  Presented  here  is  a recapitulation  of  the 
major  colloidal  forces  and  a summary  of  the  well  known  equations,  which  serves  as  a 
reference  for  colloidal  interaction  calculations  used  in  particle  deposition  modeling. 
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1-2-1  Van  der  Waals  Interactions 

The  universal  attractive  forces  between  atoms  and  molecules,  known  as  van  der 
Waals  forces,  also  operate  between  macroscopic  objects.  Van  der  Waals  attraction 
originates  from  the  instantaneous  dipole  moments  generated  by  the  temporary 
asymmetrical  distribution  of  electrons  around  atomic  nuclei. 

Essentially  there  are  two  theoretical  approaches  to  the  evaluation  of  the  van  der 
Waals  attraction.  The  microscopic  approach,  due  largely  to  Hamaker(1937),  calculates 
the  interaction  between  macroscopic  bodies  by  the  pairwise  summation  of  all  the  relevant 
intermolecular  interactions.  The  expressions  obtained  in  this  manner  may  be  split  into  a 
purely  geometric  part  and  a constant  H,  the  Hamaker  constant,  which  is  related  only  to 
the  properties  of  the  interacting  bodies  and  the  medium.  The  assumption  of  pairwise 
additivity  is  a serious  deficiency  and  always  results  an  overestimation  of  the  van  der 
Waals  force.  This  deficiency  can  be  overcome  in  the  macroscopic  approach  suggested  by 
Liftshitz  (1956),  in  which  the  interaction  is  calculated  entirely  from  considerations  of  the 
macroscopic  electromagnetic  properties  of  the  medium,  but  its  application  is  limited  to 
cases  where  detailed  knowledge  of  the  dielectric  responses  of  the  interacting  media  over  a 
wide  frequency  range  is  available.  Certain  simplifications  are  possible,  very  clearly 
described  by  Hough  and  White  (1980)  and  Prieve  and  Russel  (1988),  but  considerable 
computation  is  needed.  It  is  also  worth  noticing  that  the  difference  between  the  results 
obtained  from  the  two  approaches  is  not  too  great  ( usually  less  than  60%,  Elimelech  et 
al.  1995).  Hence,  the  practical  way  to  calculate  the  van  der  Waals  interaction  is  to  rely  on 
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the  simple  Hamaker-type  expressions  with  corrections  to  account  for  the  difference 
between  the  two  approaches  (retardation  effect).  Listed  in  Table  1-1  are  a few  expressions 
for  the  retarded  van  der  Waals  interaction  potential  between  a sphere  and  a plate,  which  is 
the  geometry  encountered  in  particle  deposition  modeling. 

1 -2-2  Double  layer  Interaction 

Practically  all  solid  surfaces  in  an  aqueous  medium  carry  a surface  charge,  which 
can  arise  from  mechanisms  such  as  ionization  of  surface  groups  and  specific  adsorption 
of  ions.  In  an  electrolyte  solution,  the  distribution  of  ions  in  the  vicinity  of  the  solid 
surface  is  determined  by  electrical  interaction  with  the  surface  and  the  thermal  motions  of 
the  ions  (diffusion).  Counter-ions  are  attracted  toward  the  surface  and  can  be  either 
closely  associated  with  the  surface  or  distributed  exponentially  into  the  solution.  Solid 
surfaces  of  the  same  charge  approaching  each  other  begin  to  experience  repulsion  when 
their  diffusion  layers  overlap.  On  the  other  hand,  solid  surfaces  with  different  charges 
undergo  an  attraction  when  they  come  close  to  each  other. 

In  considering  the  double  layer  interaction,  an  important  quantity  is  the  electric 
potential  at  the  inner  boundary  of  the  diffusion  layer,  i.e.,  close  to  the  surface  but  just 
outside  a layer  of  closely  associated  counter-ions,  usually  called  Stem  layer.  This 
potential  is  not  measurable  directly,  but  is  close  to  another  quantity,  the  zeta-potential, 
which  can  be  obtained  from  a number  of  electrokinetic  experiments  such  as  the  streaming 
potential  and  electro-osmosis  measurements.  Some  of  the  widely  used  expressions  for 
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calculating  the  double  layer  interaction  potential  for  the  sphere-plate  geometry  are  listed 
in  Table  1-2. 


Table  1-1  Retarded  van  der  Waals  interactions  for  sphere-plate  geometry 


Author 

Expression 

Validity 

Schenkel  & 
Kitchener 
(1960) 

1 

1 

X 

h<X,w/7i 

h«a 

6h  |_1  + 1 1.1 2h  / Xw 

Ho  and 
Higuchi 
(1968) 

Ha 

2.45  (O 

2.17 

aw'|2  | 059 

(K)3' 

h>A.w/7r 

h«a 

6h 

107i  \ h ) 

6071 2 

V h J 280ti3 

l h J 

Gregory 

(1981) 

-HaL532hi„(i+  ** )] 

6h  [_  Xw  k 5.32h/_ 

h«a 

Czarnecki 

(1979) 

-H 

245^ 

60re 

(h-a  > 
h2 

h+3a 

217^ 

7207T2 

^h-2a  ) 
h3 

h+4a 

059?i. 
5010 Tt3 

(h-3a  > 
h4 

h+5a 

h>Xw/7i 

l(h+2a)2> 

v(h+2a)3J 

l(h+2a)4J 

Note:  = 2nc  / a>  v , where  c is  the  velocity  of  light  and  cov  is  the  dispersion 

requency. 

A.w  is  a characteristic  wavelength  with  a value  around  lOOnm  for  most  materials. 


Table  1-2  Double  layer  interaction  potential  for  sphere-plate  geometry 


Author 

Expression 

Validity 

Hogg  et  al.  (1966) 
Wiese  & Healy 
(1970) 

Usui  (1973) 

7ras(cpf  +92) 

+ for  constan 
- for  constan 

2<P,<P2  h 

.9? +92 
t potential  a 
t charge  ap; 

fl+e"KM 

Ll — e_Kh  J 

pproxima 

?roximati 

±ln(l-e_2Kh) 

ition 

on 

Small  cpj 
h«a 
Ka>>l 

Gregory  (1975) 

nSTran^T6^  tanh  Ze(p>  tanh  ZC(P2 
B k2  4kBT  4kBT 

h«a, 

Kaj>5 

Symmetric 

electrolyte 

e Yn-  z- 

Note:  k = , the  reciprocal  Debye  length. 

V ekBT 

Zj  is  the  valence  of  ion  i. 
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1-2-3  Non-DLVO  Interactions 

Since  the  development  of  the  DLVO  theory,  many  investigators  reported 
experimental  results  which  can  not  be  explained  by  the  consideration  of  the  two  principle 
forces.  In  such  cases,  additional  interactions  have  been  proposed.  Among  them  the  most 
widely  accepted  ones  are  the  hydration,  hydrophobic,  steric  and  polymer  bridging 
interactions,  which  will  be  discussed  briefly.  It  should  also  be  noted  that  the  quantitative 
expressions  for  these  non-DLVO  forces  are  still  not  available,  therefore  they  can  not  be 
incorporated  into  the  modeling  work. 

Hydration  force 

The  nature  of  water  close  to  a particle  surface  can  be  very  different  from  bulk 
water,  for  these  portion  of  water  tend  to  be  bounded  more  tightly  to  the  surface  by  ionic 
groups  or  hydrophilic  sites  on  the  surface.  The  approach  of  a particle  toward  a collector 
with  one  or  both  of  the  surfaces  hydrated  will  generally  be  hindered  by  an  extra  repulsive 
interaction.  This  hydration  repulsion  arises  essentially  from  the  need  for  the  surfaces  to 
become  dehydrated  if  true  contact  between  the  particle  and  the  collector  surface  is  to 
occur. 

Hydrophobic  force 

When  a surface  has  no  polar  or  ionic  groups  or  hydrogen-bonding  sites,  there  is 
no  affinity  for  water  and  the  surface  is  said  to  be  hydrophobic.  The  nature  of  water  in 
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contact  with  such  a surface  will  be  different  from  that  of  the  bulk  water,  which  is 
significantly  structured  because  of  the  hydrogen  bonding  between  the  molecules.  Water 
confined  in  the  gap  between  two  such  surfaces  will  not  be  as  structured  as  the  bulk  water 
when  the  two  hydrophobic  surfaces  come  close  to  each  other.  For  a narrow  gap,  this 
would  result  in  an  increased  free  energy  of  the  water  relative  to  the  bulk  water.  Hence, 
the  water  in  the  gap  has  the  tendency  to  leave  the  gap.  In  other  words,  there  would  be  an 
attraction  between  hydrophobic  surfaces  when  they  are  close  enough. 

Steric  and  Polymer  bridging  force 

Polymer  adsorption  onto  solid  surfaces  can  greatly  change  the  interactions.  With  a 
very  small  adsorption  amount,  the  individual  polymer  chains  can  be  attached  to  two  or 
more  surfaces,  providing  an  extra  attraction  between  them.  This  effect  is  called  polymer 
bridging.  On  the  other  hand,  if  a big  amount  is  adsorbed  onto  the  surfaces,  the  polymer 
layer  can  cause  a strong  repulsive  force  between  them,  usually  termed  as  steric  repulsion. 
As  the  surfaces  approach  sufficiently  close,  the  adsorbed  layers  come  into  contact  and 
any  closer  approach  would  involve  some  penetration  of  the  hydrophilic  chains.  Since 
these  chains  are  hydrated,  overlapping  of  the  layers  would  cause  some  degree  of 
dehydration  and  hence  a repulsion  between  the  surfaces. 

1-2-4  Combined  Interaction 

The  total  interaction  energy  between  two  surfaces  is  the  summation  of  all  the 


contributions. 
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O = dhydW  + Odl  + + ^HP  + + ^PB  ■* ( 1 -26) 

Here,  Ovdw,  <bDL,  ®H,  ®HP,  Os  and  ®PB  represent  the  van  der  Waals,  double  layer, 
hydration,  hydrophobic,  steric  and  polymer  bridging  interaction  energies,  respectively. 
Up  to  now,  only  the  van  der  Waals  and  double  layer  interactions  can  be  predicted  reliably 
through  theoretical  models  as  discussed  above.  However,  with  the  development  of 
surface  interaction  force  measurement  instrumentation,  such  as  the  atomic  force 
microscopy  (AFM)  and  the  total  internal  reflection  microscopy  (TIRM),  the  total 
interaction  energy  can  be  measured  accurately.  This  information  instead  of  the 
theoretical  predictions  can  be  more  reliably  used  as  the  input  data  in  the  particle 
deposition  modeling. 

A typical  curve  for  the  combined  interaction  energy  as  a function  of  the  separation 
distance  with  consideration  of  only  the  van  der  Waals  and  double  layer  interactions  is 
shown  in  Figure  1-1.  The  most  important  features  of  this  curve  include  a deep  primary 
minimum  and  an  energy  barrier.  The  existence  of  the  deep  primary  minimum  means  that 
once  a particle  comes  into  this  region  it  will  be  dragged  toward  the  surface  and 
permanently  attached  to  it.  On  the  other  hand  particles  must  overcome  the  energy  barrier 
in  order  to  come  into  the  primary  minimum.  Because  of  the  different  distance 
dependence  of  the  van  der  Waals  and  the  double  layer  potentials,  the  former  is  always 
greater  than  the  latter  at  sufficiently  large  separation  distances.  This  gives  a secondary 
minimum  in  the  potential  curve,  which  can  be  responsible  for  particle  deposition  too. 
Since  the  interaction  potentials  are  directly  proportional  to  particle  size,  secondary 
minimum  effects  are  more  significant  with  larger  particles  (greater  than  about  1pm  in 
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diameter,  Elimelech  et  al.  1995).  The  possible  particle  deposition  in  the  secondary 
minimum  will  be  neglected  in  the  current  research  because  only  the  colloidal  particles 
will  be  treated. 


Figure  1-1  Van  der  Waals,  double  layer  and  the  combined  interaction  potentials  between 
a 0.5  pm  (diameter)  particle  and  a plate.  (The  surface  potentials  are  assumed  to  be 
16.8mV,  the  solution  is  a 20mM  1-1  electrolyte,  the  Hamaker  constant  is  8.5xl0'21  J,  the 
reciprocal  Debye  length  is  2.17nm) 


1-2-5  Effect  of  Surface  Heterogeneity  on  Surface  Interactions 

Most  of  the  work  on  surface  interactions  have  been  focused  on  ideal  surfaces: 
molecularly  smooth  and  uniform  distribution  of  properties.  In  the  recent  years,  it  becomes 
clearer  that  the  surface  heterogeneity,  such  as  non-uniform  distribution  of  surface  charge 
and  the  surface  roughness,  may  play  an  important  role  in  altering  the  particle  deposition 
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behavior.  Most  solid  surfaces  in  aqueous  media  are  heterogeneously  charged  at  both  the 
microscopic  (molecular)  and  macroscopic  levels  ( Elimelech  et  al.  1995).  Surface  charge 
heterogeneity  is  attributed  to  the  complexity  of  the  crystaline  structure  of  solids  and  to 
their  complex  chemical  composition  (Spoito  1984,  Jaroniec  and  Madey  1988).  Surface- 
bound  impurities  may  be  an  additional  source  of  charge  heterogeneity. 

Due  to  the  charge  heterogeneity  on  the  surfaces,  especially  that  on  the  collector 
surface,  some  sites  may  be  favorable  while  others  may  be  unfavorable  to  particle 
deposition  and  this  can  greatly  affect  the  overall  particle  deposition  (Hull  and  Kitchener 
1969,  Bowen  and  Epstein  1979,  Gregory  and  Wishart  1980,  Adamczyk  et  al.  1983, 
Dabros  and  van  de  Ven  1983,  Vaidyanathan  and  Tien  1991,  Kihira  et  al.  1992,  Song  et  al. 
1994). 

The  effect  of  surface  roughness  on  colloidal  interactions  has  also  been 
investigated  in  the  recent  years  through  modeling  (Elimelech  and  O’Melia,  1990a,  Suresh 
and  Waltz  1996).  Rough  surfaces  may  decrease  the  total  energy  barrier. 

1-3  Deposition  of  Colloidal  Particles  onto  Solid  Surfaces 

Deposition  of  small  solid  particles  from  a liquid  dispersion  are  encountered  in  a 
variety  of  situations  in  manufacturing  and  process  technologies  and  in  natural  aquatic 
environments.  In  some  cases,  deposition  is  desirable,  as  in  filtration  and  various 
processes  in  which  a coating  of  deposited  particles  needs  to  be  formed.  In  many  other 
instances,  deposition  needs  to  be  prevented.  These  include  detergency,  mineral 
processing,  where  slime  coating  of  colloidal  particles  onto  larger  grains  can  be  a serious 
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problem,  and  the  fouling  of  heat  exchanger  surfaces.  A number  of  biological  examples 
could  also  be  mentioned,  such  as  the  deposition  of  bacteria  on  the  tooth  surface  and 
thrombus  formation  associated  with  artificial  organs.  There  is  considerable  interest  in  the 
transport  of  colloidal  contaminants  including  bacteria  and  viruses  in  the  ground  water  and 
this  depends  greatly  on  the  deposition  behavior  of  such  particles  in  soils  and  aquifers. 
Fundamental  studies  of  particle  deposition  have  been  undertaken  with  different 
perspectives  and  with  various  objectives  (Tien  1989,  Elimelech  et  al.  1995).  The  current 
section  presents  the  key  outlines  of  the  two  general  approaches,  Lagrangian  and  Eulerian, 
used  in  the  particle  deposition  modeling,  followed  by  a detailed  description  of  the 
Eulerian  method  which  will  be  applied  in  the  current  study. 


1-3-1  Filtration  Efficiency,  Penetration,  Mass  Transfer  Coefficient  and  Filter  Coefficient 

A filter’s  ability  to  collect  particles  is  usually  expressed  in  one  of  the  following 
four  qualities:  filtration  efficiency,  penetration,  filter  coefficient  and  mass  transfer 
coefficient.  The  filtration  efficiency  r\  for  a certain  depth  L of  filter  is  defined  as 


r|  = 


cin  — cefY 


(1-27) 


where  cin  and  cefT  denote  the  influent  and  effluent  particle  concentrations,  respectively. 
The  penetration,  most  often  used  in  aerosol  science,  is  defined  as  the  effluent  and  influent 


concentration  ratio  cef1/cjn,  or  1 -rj . These  two  qualities  are  the  overall  measurement  of  the 
collection  ability.  The  filter  coefficient  A,  and  the  mass  transfer  coefficient  km  are,  on  the 
other  hand,  local  variables,  which  may  change  with  position.  But  in  most  cases  these  two 
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qualities  were  given  as  averaged  values  in  the  whole  filter  medium.  The  definition  of  the 


mass  transfer  coefficient  is 


k 


J 


(1-28) 


m 


c csurf 


where  J is  the  particle  deposition  flux,  c the  bulk  particle  concentration  and  csurf  is  the 
particle  concentration  on  the  collector  surface.  In  most  of  the  particle  deposition  studies, 
a perfect  sink  of  particles  is  assumed  on  the  collector  surface,  i.e.,  csurf=0.  That  is  because 
the  particles  themselves  become  part  of  the  collector  once  they  are  attached  to  the 
collector  surface.  Balancing  the  movement  of  particles  in  a differential  depth  of  filter 
gives 


filter,  1 is  the  coordinate  of  the  filter  depth.  If  k,„  is  assumed  to  be  constant  along  the  depth 
of  the  filter,  substituting  equation  (1-28)  into  (1-29)  and  carrying  out  the  integration 
yields 


which  shows  that  the  particle  concentration  decays  exponentially  in  the  direction  of  the 
flow.  The  exponential  constant  is  called  filter  coefficient,  X,  i.e., 


(1-29) 


where  U is  the  superficial  velocity  of  the  suspension,  Asp  the  specific  surface  area  of  the 


c = c0e  u 


(1-30) 


(1-31) 


The  relationship  between  the  filtration  efficiency  and  the  filter  coefficient  is 
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r|  = 1 - e_^L  (1-32) 

It  is  obvious  that  the  mass  transfer  coefficient  is  the  most  intrinsic  property  of  the  particle 
deposition  process,  therefore  km  will  be  used  throughout  the  current  work. 


(a)  Lagrangian  Approach 


(b)  Eulerian  Approach 


Figure  1-2  Two  approaches  for  particle  deposition  modeling 


1-3-2  Two  Approaches  in  Particle  Deposition  Modeling 

Due  to  its  practical  importance,  numerous  studies  have  been  conducted  to  model 
the  particle  deposition  process  during  the  past  several  decades  (Rajagopalan  and  Tien 
1979,  Tien  1989,  Elimelech  et  al.  1995).  Two  principal  approaches  are  available  for 
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modeling  purposes.  The  Lagarangian  approach  considers  the  behavior  of  individual 
particles,  or  their  trajectories,  at  the  microscopic  level,  whereas  the  Eulerian  approach 
describes  the  particles  collectively  in  terms  of  their  distribution,  or  probability  density,  in 
time  and  space,  thus  it  is  a typical  macroscopic  approach.  The  general  procedures  of 
these  two  basic  methodologies  are  summarized  in  Figure  1-2.  In  the  current  research,  we 
will  take  the  Eulerian  approach,  because  all  the  equations  derived  by  using  this  method 
can  be  directly  applied  to  other  heat  and  mass  transfer  processes  in  the  same  medium  due 
to  the  similarity  between  them. 

1-3-3  General  Procedure  for  the  Eulerian  Approach 

Figure  1-3  is  an  exploded  detail  of  the  Eulerian  method.  Also  shown  in  the  same 
figure  are  the  two  major  tasks  of  the  current  research:  an  improved  porous  medium 
model  and  the  particle  deposition  modeling  based  on  it. 

Removal  of  small  colloidal  particles  in  a porous  filter  is  not  achieved  by  straining 
since  the  pore  size  of  the  porous  filter  is  usually  a few  order  of  magnitude  larger  than  the 
colloidal  particles  being  removed.  Thus  the  only  way  for  the  particles  to  be  removed  is 
by  deposition  on  the  collector  surfaces  through  colloidal  interactions  between  particles 
and  collectors.  Deposition  of  particles  from  an  aqueous  fluid  medium  flowing  through  a 
porous  medium  may  be  viewed  as  a two-step  process:  the  transport  of  the  particles  from 
liquid  bulk  to  the  proximity  of  the  fiber  surfaces  by  the  combination  of  convection  and 
Brownian  diffusion  (which  we  may  call  “bulk  transport”)  and  the  subsequent  adhesion  of 
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Figure  1-3  General  procedure  for  modeling  of  particle  deposition  rate  in  fibrous  or 
granular  media  used  in  the  current  research 


particles  to  the  collector  surface  by  colloidal  forces  (“surface  deposition”).  The 
interaction  forces  of  non-hydrodynamic  origin  include  van  der  Waals  attraction  and 
electrostatic  (or  double  layer)  interactions  for  most  particulate  systems.  In  some  cases, 
hydrophobic,  steric,  and  hydration  forces  also  play  an  important  role  especially  for 
biological  particles  or  the  particles  with  absorbed  polymers  or  surfactants.  Since  all  these 
surface  forces  are  very  short  ranged,  they  usually  take  effect  only  when  the  particles  come 
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very  close  to  the  collector  surface  (<  0.1  (am).  Since  the  colloidal  particles  are  much 
smaller  than  the  collector,  the  interactions  between  particles  and  the  collectors  may  be 
considered  as  interactions  between  particles  and  a flat  surface. 

The  concentration  distribution  of  particles  convected  with  a fluid  passing  through 
the  porous  medium  can  be  described  by  the  following  steady  state  convective  diffusion 
equation.  It  includes  the  contribution  of  an  external  force  field  (Bird  et  al.  1959; 
Spielman  and  Friedlander  1 974) : 


u • Vc  = V • 


£>Vc  + 


V 

kaT 


VOc 


(1-33) 


Here  u is  the  fluid  velocity,  c the  particle  concentration,  *D  the  particle  diffusion 
coefficient,  O the  total  interaction  potential,  kB  the  Boltzmann  constant  and  T the  absolute 
temperature.  Since  the  colloidal  particles  are  so  small,  their  size  is  neglected  and  the 
particle  laid  fluid  is  viewed  as  a continuum.  When  the  particles  are  very  close  to  the 
collector  surface,  however,  the  size  effect  may  become  significant  as  the  diffusivity  is 
influenced  significantly.  This  hydrodynamic  interaction  at  close  range  is  usually  taken 
care  of  by  lumping  its  contribution  into  a position-dependent  diffusion  coefficient: 
£>(h)  = 2 / g(h) . Here  is  the  particle  diffusion  coefficient  in  the  bulk,  which  may  be 
estimated  by  the  Stokes-Einstein  equation  (Bird  et  al.  1960).  g(h)  is  the  hydrodynamic 
retardation  function  g(h)  given  by  Brenner  (1961). 

Although  colloidal  particles  are  small,  their  diffusivity  is  still  very  small 


compared  to  those  of  molecular  species.  Therefore,  convective  transport  is  usually 
dominant  over  diffusion  (i.e.,  Peclet  number  Pe»  1)  even  when  the  Reynolds  number  is 
small.  Consequently,  a diffusion  boundary  layer  exists.  In  addition,  in  many  cases  of 
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practical  interest,  the  thickness  of  the  electrical  double  layer  (which  characterizes  the 
range  of  electrostatic  interaction)  is  even  smaller  than  the  diffusion  boundary  layer 
thickness.  In  this  case,  the  prescribed  problem  can  be  simplified  in  which  the  external 
force  term  in  equation  (1-33)  can  be  dropped  while  the  contribution  of  the  colloidal 
interactions  are  represented  by  a modified  boundary  condition  on  the  collector  surface: 

= 0-34) 

Here  n is  the  coordinate  in  the  outward  normal  direction  from  the  collector  surface,  and 
K0  is  a constant  which  depends  on  the  particle-surface  interaction  potential  as  follows: 

=Z’.{fo”[s(h)«®'kBT-l]tiy}’1  0-35) 

where  h is  the  separation  distance  between  a particle  and  the  surface  and  h=y-a.  This 
simplification  is  often  called  the  Interaction  Force  Boundary  Layer  (IFBL)  approximation 
which  was  first  described  independently  by  Ruckenstein  & Prieve  (1973)  and  Spielman 
& Friedlander  (1974).  Later,  this  IFBL  approximation  was  rigorously  justified  by  the 
method  of  matched  asymptotic  expansions  (Shapiro  et  al.  1990). 

As  we  may  note  from  equation  (1-33),  the  flow  field  u needs  to  be  specified  in 
order  to  determine  the  particle  deposition  rate  on  to  the  collector  surface.  As  it  was 
pointed  out,  the  diffusive  transport  and  subsequent  deposition  of  particles  occurs  within 
the  diffusion  boundary  layer.  Thus,  the  flow  field  description  is  required  only  in  the 
immediate  vicinity  of  the  collector  surface  in  solving  equation(l-33).  Therefore,  the 
linearized  velocity  field  as  in  typical  boundary  layer  equations  will  suffice  the  need.  For 
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the  cylindrical  and  spherical  geometries,  the  linearized  velocity  fields  are  represented  by 
the  following  stream  functions,  respectively: 


\\i  = -2Ac(r  - l)2  sin0 

(1-36) 

v| f = -^  As(r  - 1)2  sin2  0 

(1-37) 

Here  Ac  and  As  are  the  flow  field  parameters  which  are  unique  to  each  of  the  porous 

medium  models. 

Using  the  IFBL  approximation,  Spielman  and  Friedlander  (1974)  obtained  the 

following  particle  deposition  rates  on  a cylindrical  collector  and 

a spherical  collector 

using  the  linearized  flow  fields: 

Shc=0.731Ai'3Pe'/3Tc(Pc) 

(1-38) 

Shs=0.635AfPe‘'3ys(Ps) 

(1-39) 

where 

pc  = 1.022A'l/3Pe‘l/3Da 

(1-40) 

ps  = 1.125Ay3Pe-,/3Da 

(1-41) 

Here  Sh=  l^a IV,  Pe=  Ua IV  and  Da=  K^a IV.  V is  the  particle  diffusivity,  a the  collector 
radius,  K^,  the  virtual  reaction  constant.  yc(Pc)  and  ys(Ps)  are  the  functions  tabulated  in 
Table  1-3.  If  the  virtual  surface  reaction  is  instantaneous,  y (P) =1  in  equations  (1-38) 
and(l-39),  the  equations  are  then  reduced  to  the  ones  without  reaction  on  the  surface, 
which  were  first  derived  by  Levich  (1962).  For  an  isolated  sphere,  the  flow  field  is 
given  by  the  Stokes  equation,  As=l.  So,  the  flow  field  parameter  As  in  equation  (1-38) 
represents  the  influence  of  neighboring  spheres  in  a granular  packing. 
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Table  1-3  Numerical  values  of  the  functions  yc  and  ys 
(converted  from  Spielman  and  Friedlander,  1974) 


Pc  or  Ps 

Yc(Pc) 

Ys(Ps) 

0 

0 

0 

0.01 

0.014044 

0.013198 

0.02 

0.02772 

0.026078 

0.05 

0.066686 

0.062857 

0.1 

0.125445 

0.118709 

0.2 

0.2241 

0.213217 

0.5 

0.423033 

0.407133 

1 

0.5984 

0.5822 

2 

0.7518 

0.738867 

5 

0.885083 

0.877833 

10 

0.939455 

0.935364 

20 

0.968952 

0.966667 

50 

0.987353 

0.986373 

100 

0.993663 

0.993069 

00 

1 

1 

1-4  Interaction  Force  Boundary  Layer  Approximation 

As  shown  above,  the  particle  deposition  rate  can  be  found  by  solving  the 
convective  diffusion  equation  (1-33)  under  the  surface  interaction  force  field  if  the  flow 
field  is  known.  But  this  equation  is  difficult  to  solve  analytically.  An  approximate 
analytical  solution  for  calculating  the  particle  deposition  rate  when  electrical  double  layer 
repulsive  forces  are  significant  may  be  useful.  The  most  widely  used  one  is  the  so  called 
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Interaction  Force  Boundary  Layer  (IFBL)  approximation.  The  IFBL  approximation  was 
developed  by  Ruckenstein  and  Prieve  (1973)  and  Spielman  and  Friedlander  (1974) 
independently.  The  basic  ideas  behind  the  approximation  are  presented  here. 

Colloidal  interaction  forces  take  effect  only  over  very  short  distances  from  the 
surface  of  the  collector,  usually  less  than  0.1  pm.  It  is  possible  to  identify  two  regions 
adjacent  to  the  collector  surface,  characterized  by  the  relative  magnitude  of  the  colloidal 
interaction  forces.  In  the  vicinity  of  the  collector  surface,  the  colloidal  interaction  forces 
predominate  while  the  convective  transport  may  be  neglected.  It  is  this  region  that  is 
referred  to  as  the  interaction  force  boundary  layer.  The  thickness  of  this  region  is 
characterized  by  the  operating  distance  of  the  double  layer,  on  the  order  of  the  Debye 
length  k'1.  In  the  region  away  from  the  IFBL,  the  colloidal  interaction  forces  vanish 
while  convective  transport  becomes  important.  The  thickness  of  this  region  is  on  the 
order  of  the  diffusion  boundary  layer,  SD,  of  the  order  aPe'l/3(Levich  1962). 

The  IFBL  approximation  is  based  on  the  assumption  that  the  diffusion  boundary 
layer  is  very  thick  compared  to  the  IFBL,  i.e., 


In  the  IFBL  region,  where  the  fluid  velocity  is  negligible,  the  transport  of  particles  is 
controlled  by  diffusion  and  colloidal  forces.  The  convective  diffusion  equation  (1-32)  can 
be  simplified  to: 


KaPe“1/3  » 1 


(1-42) 


cV 

V-  £>Vc  + — — VO  =0 
V koT  ) 


(1-43) 


With  the  assumptions  that  the  total  interaction  potential  d>  is  only  a function  of  the 
boundary  layer  coordinate  y,  which  is  perpendicular  to  the  collector  surface,  and  that  the 
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tangential  diffusion  is  negligible  compared  to  the  diffusion  perpendicular  to  the  surface, 
equation  (1-43)  can  be  easily  solved  to  get: 


J(x)  * 
c = — — e 


<t> 

' o' 

kBT 

fy 

JO 

g(h)ekBT-l 

dy  + y 

(1-44) 


This  equation  describes  the  particle  concentration  profile  in  the  IFBL  as  a function  of  the 
distance  away  from  the  surface,  which  is  shown  schematically  in  Figure  1-4. 


Figure  1-4  Concentration  profile  in  the  vicinity  of  the  collector  surface 


At  distances  much  larger  than  the  thickness  of  the  double  layer  (y»K'')  but  still 
close  enough  to  the  collector  surface,  both  the  convection  and  the  surface  interaction 
forces  may  be  neglected.  In  this  instance,  equation  (1-43)  can  be  simplified  to: 


c = — 


J(x) 


( o > 

f 00 

Jo 

g(h)ek«T-l 

dy  + y 

< ) 

(1-45) 
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which  is  linear  in  y as  shown  in  Figure  1-4.  In  this  equation,  the  upper  limit  of  the 


integral  has  been  arbitrarily  chosen  to  be  oo,  because  J'£°^g(h)e<I>/kBT  -ljdy  = 0 for  any 

value  b»k‘‘,  and  J^g(h)e°/kBl  - ljdy  becomes  a constant. 

At  distances  further  away  from  the  collector  surface,  equation  (1-33)  is  reduced 


to: 


u • Vc  = £>  V2c 


(1-46) 


It  is  obvious  that  the  linear  portion  of  the  solution  (1-44),  i.e.,  equation  (1-45)  is  the 
solution  of  (1-46).  Taking  any  point  along  this  linear  part  as  one  of  the  boundary 
conditions  to  (1-46),  we  can  then  solve  it  to  obtain  the  concentration  profile  outside  of  the 
IFBL.  The  solution  we  get  would  be  the  same  if  we  use  the  extrapolated  point  at  y=0, 


♦ 

c = 


J(x) 


r 

Jo 


cD  ' 

g(h)ekBT-l 


dy 


(1-47) 


as  the  boundary  condition  instead  of  an  arbitrary  point  as  described  above. 

Equation  (1-47)  shows  that  the  problem  of  deposition  of  colloidal  particles  with  a 
potential  barrier  confined  to  a narrow  region  close  to  the  collector  surface  can  be  found 
by  solving  the  usual  convective  diffusion  equation  (without  the  external  force  field) 
subject  to  the  following  boundary  condition: 


J(x)  = -K0c  (x) 


(1-48) 


where 


K*  = 


Dn 


J”(g(h)e0/kBT)dy 

The  boundary  condition  (1-48)  is  equivalent  to  a first  order  reaction  on  the  surface. 


(1-49) 
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1-5  Dissertation  Organization 

As  stated  in  Section  1-3,  the  major  tasks  of  the  current  research  are  two  folds:  an 
improved  porous  medium  model,  the  effective  medium  approximation,  which  will  be 
applied  to  both  fibrous  and  granular  media  and  will  also  be  extended  to  packings  of 
polydisperse  cylinders  and  spheres;  and  the  particle  deposition  modeling  based  on  the 
improved  porous  medium  model. 

Chapter  1 provides  the  necessary  fundamentals  for  understanding  the  topic  and 
also  the  major  research  developments  in  this  field  which  will  be  frequently  referred  to  in 
the  dissertation.  This  includes  a brief  summary  of  the  porous  medium  models  and  the 
corresponding  flow  field  and  permeability  predictions,  the  surface  interaction  forces  and 
the  corresponding  equations  to  calculate  their  potentials,  and  an  outline  of  the  procedures 
for  particle  deposition  modeling,  followed  by  a brief  derivation  of  the  interaction  force 
boundary  layer  approximation  which  is  the  key  simplification  toward  solving  the 
complicated  convective  diffusion  problem. 

For  simple  regular  array  of  cylindrical  cylinders,  the  numerical  simulation  can 
provide  accurate  predictions  of  the  detailed  flow  and  concentration  fields  by  solving  the 
Stokes  and  convective  diffusion  equations.  This  is  done  in  Chapter  2. 

The  analyses  of  the  effective  medium  approximation  applied  to  randomly  packed 
fibrous  and  granular  media  are  elaborated  in  Chapter  3 and  Chapter  4,  respectively.  The 
Stokes  and  Brinkman’s  equations  are  solved  to  get  the  representative  velocity  fields.  The 
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packing  permeabilities  are  then  obtained  through  force  balances  in  the  media.  The  results 
are  compared  with  the  predictions  by  other  models  and  experimental  data. 

Chapter  5 treats  the  particle  deposition  in  fibrous  and  granular  media  based  on  the 
flow  fields  obtained  in  the  preceding  chapters.  Chapter  6 is  an  investigation  of  the 
influencing  parameters  such  as  electrolyte  concentration,  surface  potentials  and 
hydrodynamic  interaction  on  the  particle  deposition  rate. 

It  should  also  be  noted  that  although  the  current  research  is  focused  on  the 
deposition  of  colloidal  particles  in  fibrous  and  granular  media,  the  concluding  equations 
can  be  applied  to  other  mass  and  heat  transfer  processes  in  the  same  media. 


CHAPTER  2 

DEPOSITION  OF  COLLOIDAL  PARTICLES  IN  PERIODIC 
ARRAYS  OF  CYLINDERS 

2-1  Periodic  Array  of  Cylindrical  Collectors 

As  discussed  in  Chapter  1 , the  velocity  field  around  each  of  the  collectors  has  to 
be  obtained  in  order  to  determine  the  deposition  rate  in  a porous  medium  by  solving  the 
convective  diffusion  equation.  For  simpler  geometries,  like  periodic  arrays  of  cylindrical 
collectors  which  will  be  treated  in  this  Chapter,  some  numerical  algorithms  can  be  used  to 
solve  the  Stokes  equation  and  the  convective  diffusion  equation  simultaneously.  The 
Finite  Element  Method  (FEM)  will  be  used  in  the  current  study. 

Figure  2-1  describes  a square  array  of  collectors  which  consists  of  two  different 
size  cylinders  of  r,  and  r2  alternating  with  the  minimum  surface-to-surface  distance  of  s. 
While  all  three  collector  parameters  can  be  varied,  there  remain  two  degrees  of  freedom  if 
the  volume  fraction  of  collectors  (or  porosity)  is  to  be  fixed.  Thus  for  a fixed  value  of 
porosity,  the  ratio  of  r,/r2  can  be  varied  and  its  influence  on  the  filtration  efficiency  can  be 
investigated.  In  the  conventional  modeling  where  the  cell  models  of  Happel  or  Kuwabara 
are  used,  the  effect  of  distribution  of  collector  size  on  filtration  efficiency  cannot  be 
investigated  as  they  take  into  account  only  the  representative  size  and  the  volume  fraction 
of  collectors  (Speilman  1977,  Russel  et  al.  1989,  Tien  1989). 


31 


32 


Throughout  the  calculation,  r,  was  fixed  at  0.5  mm  while  r,  was  varied  between 
0.05  and  0.5  mm.  s was  also  varied  so  that  the  porosity  would  be  in  the  range  of  50%  to 
70%.  All  other  parameters  were  specified  to  simulate  the  filtration  of  silica  or  styrene 
latex  particles  which  are  smaller  than  1.5  pm  in  diameter  and  dispersed  in  an  aqueous 
solution.  Thus,  the  Hamaker  constant  was  set  between  10'21  and  10'19  J.  Assuming  a 
spherical  shape,  the  diffusivity  of  the  particles  was  calculated  using  the  Stokes-Einstein 
equation.  The  surface  potential  may  vary  depending  on  the  hydronium  ion  concentration 
in  the  solution,  and  a broad  range  of  values  was  used  up  to  ±60  mV  for  the  calculation. 
While  the  surface  potentials  of  collectors  and  particles  may  vary  independently,  the  same 
value  was  assigned  to  both  collectors  and  particles  for  simplicity.  This,  however,  does 
not  restrict  any  conclusions  drawn  in  the  present  study. 

Also  described  in  Figure  2-1  is  the  flow  domain  for  the  present  numerical 
calculation  which  consists  of  four  unit  cells  in  series.  The  flow  direction  is  from  left  to 
right,  and  the  equations  (1-2,  1-46)  subject  to  the  boundary  condition  (1-34)  are  solved  by 
a finite  element  method  to  determine  the  velocity  and  the  concentration  fields.  The 
particle  concentration  at  the  inlet  of  the  flow  domain  is  specified  as  1.0  uniformly  across 
the  array  as  the  concentration  is  scaled  by  the  inlet  concentration.  The  velocity  profile  at 
the  inlet  is  also  assumed  to  be  uniform  since  it  is  not  known  a priori.  This  flat  velocity 
profile  relaxes  very  quickly  to  a new  steady  value  within  the  first  unit  cell,  and  the 
velocity  field  in  the  subsequent  three  cells  of  the  down  stream  matches  one  another 


exactly. 
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flow 


calculation  domain 


Figure  2-1  Periodic  array  of  cylindrical  collectors  and  the  flow  domain  for  numerical 
calculation 


Figure  2-2  Finite  element  mesh  used  in  the  current  investigation 
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Figure  2-3  The  velocity  field  calculated  by  FEM — stream  line 


Figure  2-4  The  particle  concentration  distribution  calculated  by  FEM — concentration 
contour  lines  (U=0.0001mm/s,  K^lxlO'9  m/s,  the  concentration  from  the  middle  toward 
outer  boundary  c/co=0. 769, 0.8, 0.83, 0.86, 0.89,0. 92, 0.95, 0.98) 
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Shown  in  Figure  2-2  2-3  and  2-4  are  the  finite  element  mesh,  the  velocity  field 
and  the  particle  concentration  distribution  calculated  by  the  Finite  Element  Method.  The 
unit  cell  described  in  Figure  2-2  contains  500  quadrilateral  elements  with  their  sizes 
decreasing  gradually  toward  the  collector  surface.  Thus,  the  flow  domain  of  calculation 
contains  a total  of  2,000  elements.  The  accuracy  of  the  present  calculation  has  been 
checked  by  comparing  the  drag  force  imposed  on  a collector  in  the  square  array  with  the 
results  of  Sangani  & Acrivos  (1982).  They  have  determined  the  drag  on  a cylinder  in  a 
square  array  by  a boundary  collocation  method  using  an  infinite  series  solution  for 
vorticity  and  stream  function  for  a low  Reynolds  number  flow.  The  dimensionless  drag 
force  on  a unit  length  cylinder  (F/pU)  was  calculated  to  be  102.90  and  532.55  when  the 
porosity  was  70%  and  50%,  respectively.  These  values  are  exactly  same  as  those  of 
Sangani  & Acrivos  to  the  second  decimal  points.  Although  there  are  no  data  available  for 
comparison,  the  calculated  concentration  field  is  expected  to  be  as  accurate  as  the 
velocity  field.  It  should  be  noted  that  the  superficial  velocity  was  set  to  an  extremely  low 
value  in  order  to  show  the  particle  concentration  change  pattern  in  the  flow  path. 

In  all  calculations  the  superficial  velocity  of  the  fluid  has  been  limited  up  to  0.1 
mm/s.  This  value  is  rather  small  compared  to  those  in  actual  filtration  processes. 
Nevertheless,  this  limit  was  necessary  in  order  to  obtain  a detectable  concentration 
change  over  the  calculation  domain  (i.e.,  four  unit  cells).  If  the  superficial  velocity  is 
higher,  the  concentration  change  over  a unit  cell  is  very  small  requiring  a longer 
calculation  domain  which  is  numerically  too  intensive.  Even  at  a velocity  as  low  as  0.001 
mm/s,  however,  the  Peclet  number  (Pe  = Ua I'D)  is  still  as  large  as  10\  and  consequently 
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the  convection  is  still  dominant  in  accordance  with  the  assumptions  described  in  Section 
1-4.  Thus,  all  physical  phenomena  occurring  at  a fluid  velocity  of  order  1 mm/s  should 
be  preserved  even  at  a fluid  velocity  as  low  as  0.001  mm/s.  Unless  the  fluid  velocity  is 
unrealistically  large  (e.g.,  1 m/s),  the  Reynolds  number  is  always  smaller  than  about  10'3 
justifying  the  use  of  Stokes  equation  for  the  velocity  field. 


2-2  Effect  of  Surface  Potentials 

In  Figure  2-5,  the  mass  transfer  coefficient  is  plotted  against  the  surface  potentials 
under  the  conditions  specified  in  the  figure  caption.  When  the  product  of  particle  and 
collector  surface  potentials  is  higher  than  about  130  (mV)2,  the  mass  transfer  coefficient 
is  very  small  no  matter  whether  the  fluid  velocity  is  high  or  low.  It  is  due  to  the  strong 
repulsion  (or  a high  repulsive  barrier)  between  the  particles  and  the  collectors  which 
prevents  the  particle  capture.  As  the  surface  potential  is  decreased,  the  repulsive  barrier 
is  lowered  and  the  number  of  particles  carried  over  the  barrier  increases  resulting  in  a 
gradual  increase  in  mass  transfer  coefficient.  When  the  product  of  the  surface  potentials 
is  lower  than  about  100  (mV)2,  the  rate  of  particle  capture  does  not  increase  any  more 
since  it  is  now  limited  by  the  diffusion  which  brings  the  particles  to  the  proximity  of 
collector  surfaces.  It  may  be  interesting  to  note  that  the  critical  value  of  the  product  of 
the  surface  potentials  at  which  the  particle  capture  is  diffusion  limited  is  positive.  It 
implies  that  collectors  with  opposite  charge  to  the  particles  may  not  improve  the  filtration 
efficiency  any  further  since  the  particle  capture  is  still  diffusion  limited. 
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Figure  2-5  Dependence  of  mass  transfer  coefficient  on  surface  potentials.  (0.5  pm 
particles  in  diameter;  temperature  = 298K;  Hamaker  constant  = 8.5  x 10'2IJ;  ionic 
concentration  = 2mM,  1-1  electrolyte;  r,=0.25mm,  r,-0.5  mm,  s=0.25mm). 


The  existence  of  a small  region  of  the  critical  surface  potential  product  below 
which  the  particle  deposition  is  limited  by  diffusion  and  beyond  which  it  is  essentially 
impossible  was,  in  fact,  predicted  by  previous  investigators  (Ruckenstein  and  Prieve 
1973)  and  also  observed  experimentally  (Spielman  and  FitzPatrick,  1973).  This  will  be 


explored  in  detail  later  in  Chapter  6. 


38 


Figure  2-6  Influence  of  collector  size  ratio  and  porosity  (U  = 0.01  mm/s;  1^=0.001  m/s; 
D=8. 73x1 0-1 3m2/s) 


2-3  Effect  of  Collector  Size  Ratio  and  Porosity 

In  Figure  2-6,  the  influence  of  porosity  and  the  r,/r2  ratio  on  the  mass  transfer 
coefficient  is  shown.  As  we  may  anticipate,  the  mass  transfer  coefficient  decreases  with 
increasing  porosity.  It  is  due  to  the  more  favorable  flow  (or  streamline)  pattern  for 
particle  capture  created  at  a smaller  r,/r2  ratio.  When  r,/r2  = 1 (see  Figure  2-1),  the 
collector  assemblage  is  simply  a square  array  of  same  size  cylinders.  When  r,  is 
negligibly  small  compared  to  r2,  on  the  other  hand,  the  assemblage  is  virtually  a staggered 
array  which  should  increase  the  capture  efficiency  as  the  flow  pattern  in  this  array 
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decreases  the  thickness  of  the  mass  transfer  boundary  layer.  We  may  anticipate  that  this 
dependence  of  mass  transfer  coefficient  on  r,/r2  ratio  becomes  more  significant  with 
decreasing  porosity  as  predicted  by  the  present  calculation. 


0 0.2  0.4  0.6  0.8  1 1.2  1.4  1.6 

Particle  Diameter  (Pm) 


Figure  2-7  Dependence  of  mass  transfer  coefficient  on  particle  size  (Temperature  = 

298K;  Hamaker  constant  =8.5  xlO'21  J;  ionic  concentration  = 2mM  of  1-1  electrolyte; 
r,=0.25mm,  r2=0.5mm,  s=0.25mm;  cp,=(p2=-10mV). 

The  dependence  of  mass  transfer  coefficient  on  the  array  (or  packing)  structure  of 
collectors  is  rather  significant  as  shown  in  Figure  2-6.  Such  dependence,  however, 
cannot  be  predicted  by  the  conventional  isolated-collector  models  utilizing  the  cell 
models  of  either  Happel  or  Kuwabara  since  only  a representative  collector  size  and  the 
porosity  are  taken  into  account  in  these  models  (Spielman  1977). 
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2-4  Effect  of  Particle  Size 

The  particle  size  has  significant  influence  on  the  mass  transfer  coefficient  since 
the  filtration  process  is  typically  diffusion  limited.  In  Figure  2-7  the  mass  transfer 
coefficient  is  plotted  as  a function  of  particle  size  for  which  the  diffusion  coefficient  has 
been  calculated  using  the  Stokes-Einstein  equation  for  spherical  particles.  With 
increasing  particle  size  the  mass  transfer  coefficient  decreases  due  to  smaller  diffusivity. 
When  the  particle  size  is  about  1.2  pm,  the  mass  transfer  coefficient  is  virtually  zero  since 
the  particle  diffusivity  is  very  small.  As  the  particle  size  is  further  increased,  the  mass 
transfer  coefficient  is  expected  to  increase  gradually  since  the  inertia  effect,  which  has  not 
been  considered  in  the  present  study,  will  start  to  play  a role.  Consequently,  particles 
within  the  size  range  of  about  1 to  5 pm  are  likely  to  be  most  difficult  to  capture  since 
neither  Brownian  nor  the  inertia  effect  is  significant. 


2-5  Summary 

We  have  studied  the  influence  of  the  packing  structure,  surface  potentials  of  the 
filter  and  particles,  and  particle  size  on  the  capture  efficiency  of  colloidal  particles  in  an 
idealized  square  array  of  cylindrical  collectors  using  a finite  element  method.  The 
collector  array  which  consists  of  cylinders  of  two  different  sizes  has  been  considered  for 
various  values  of  porosity.  For  fibrous  media  which  are  composed  of  simple  regular 
array  of  fibers,  the  numerical  solution  of  the  Stokes  equation  and  the  convective  diffusion 
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equation  gives  very  accurate  predictions  of  the  flow  field  and  concentration  distributions. 
The  numerical  calculation  is  a very  important  technique  for  heat  and  mass  transfer 
processes  in  porous  media  because  the  model  approach  as  the  one  discussed  in  the  next 
two  chapters  can  not  distinguish  the  detailed  structures  (different  regular  arrays)  and  the 
orientations  of  the  structure  relative  to  the  superficial  velocity.  The  calculated  results  can 
not  only  be  used  for  predictions  of  deposition  of  colloidal  particles  in  ideal  fibrous  media 
but,  more  importantly,  can  be  extended  to  predict  the  heat  and  mass  transfer  in  laminar 
flow  though  fibrous  and  granular  media  due  to  the  analogy  between  these  processes. 
With  the  development  of  computer  technology  and  the  software  for  numerical  analysis, 
more  and  more  complicated  structures  will  be  able  to  be  attacked  efficiently. 


CHAPTER  3 

EFFECTIVE  MEDIUM  APPROXIMATION  FOR  FIBROUS  MEDIA 
3-1  Concept  of  Effective  Medium  Approximation 

For  many  applications,  such  as  in  modeling  of  particle  removal  (or  filtration)  and 
heat  transfer  in  porous  media,  a good  model  should  describe  well  not  only  the 
macroscopic  behavior  of  the  medium  but  also  the  flow  characteristics  around  the 
individual  element  (or  collector  in  case  of  particle  removal  process)  in  a smaller  length 
scale  because  the  flow  field  around  the  collector  surface  affects  these  processes 
significantly.  Both  the  cell  models  and  the  Brinkman’s  model  provide  a flow  field  around 
a representative  packing  element,  and  have  been  used  to  predict  the  particle  removal  rate 
in  the  fibrous  and  granular  media  (Tien  1989).  While  the  cell  models  appear  to  provide 
good  predictions  for  permeability,  the  Brinkman’s  model  tends  to  underestimate  the 
permeability  when  the  porosity  is  low.  However  it  should  be  pointed  out  that  a good 
prediction  of  the  permeability  of  a medium  is  a prerequisite,  and  doesn’t  necessarily 
results  in  a good  prediction  of  filtration  rate  or  heat  transfer  rate.  In  fact,  it  was  realized 
recently  (Wang  and  Sangani,  1997)  that  very  different  flow  fields,  which  result  in 
different  heat  and  mass  transfer  rate  in  a fibrous  medium,  may  yield  the  same 
permeability.  So  the  validity  of  the  flow  fields  given  by  these  models  is  still  an  open 
question. 
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Effective 


Fluid 

Figure  3-1  Schematic  of  a model  unit  for  the  Effective  Medium  Approximation  (EMA) 

In  the  present  work,  another  model  for  porous  media  is  presented  in  which  the 
effective  medium  approximation  (EMA)  is  applied.  This  model  is,  in  a sense,  a 
combination  of  the  cell  model  and  Brinkman’s  model  as  it  assumes  a representative 
packing  element  surrounded  by  a fluid  envelope  up  to  a certain  radius  and  by  an  effective 
medium  beyond  it  (Figure3-1).  While  this  model,  like  the  Brinkman’s  model,  may 
account  for  the  influence  of  neighboring  elements  on  the  flow  better  than  cell  models,  the 
tendency  of  Brinkman’s  model  to  underestimate  the  permeability  is  apparently  corrected 
by  the  presence  of  a fluid  envelope  in  this  model.  In  addition,  as  it  will  be  shown,  this 
model  provides  more  accurate  prediction  of  particle  removal  rate  (or  heat  transfer  rate) 
compared  to  other  models. 

The  effective  medium  approximation  was  successfully  used  recently  by  Dodd  et 
al.  (1995)  in  the  study  of  hydrodynamic  interactions  on  diffusivities  of  integral  membrane 
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proteins.  It  has  been  also  applied  in  estimating  the  permeability  of  flow  through  an 
parallel  array  of  randomly  packed  cylinders  (Wang  and  Sangani  1997). 

The  effective  medium  approximation  is  essentially  a combination  of  the 
Brinkman’s  model  and  the  cell  model.  In  this  model,  a representative  packing  element 
(which  we  will  also  call  collector  in  particle  removal  process)  is  assumed  to  be 
surrounded  by  a liquid  envelope  up  to  a certain  radius  and  by  an  effective  medium 
beyond  it  (see  Figure  3-1).  The  liquid  envelope  is  intended  to  represent  the  environment 
around  the  collector  in  a length  scale  equivalent  to  collector  size,  whereas  the  effective 
medium  away  from  the  collector  represents  the  environment  in  a macroscopic  scale, 
thereby  accounting  for  influence  of  the  neighboring  collectors  on  the  flow  field.  In  this 
model,  the  flow  field  is  described  by  Stokes  equation  within  the  fluid  envelope  and  by  the 
Brinkman  equation  in  the  region  of  effective  medium: 

Vp=pV2u  a<r<aR  (3-1) 

Vp  = pV2u --jpu  r>aR  (3-2) 

Here,  k is  the  Darcy’s  permeability,  p fluid  viscosity,  r the  radial  coordinate  with  its 
center  at  the  origin  of  the  representative  collector,  and  a and  aR  are  the  radii  of  the 
packing  element  and  the  fluid  envelope,  respectively.  One  of  the  most  important 
parameters  in  the  EMA  model  may  be  the  size  of  the  fluid  envelope,  aR.  While  it  may 
depend  on  the  detailed  structure  of  the  cylindrical  collector  layout  around  the 
representative  collector,  the  simplest  choice  may  be  R=l/(l-e)l/n,  where  e is  the  porosity 
of  the  porous  medium  and  Q the  dimension  parameter  of  the  flow  field,  Q-2  for  a flow 
through  fibrous  media  and  Q=3  for  a flow  through  granular  media.  This  choice  of  R is 
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equivalent  to  matching  the  void  fraction  of  the  model  system  with  the  porosity  of  the 
macroscopic  medium.  We  may  note  that  the  Brinkman’s  model  (Brinkman,  1947)  can  be 
considered  as  a special  case  of  the  EM  A,  with  R=l,  which  is  equivalent  to  assuming  that 
the  representative  particle  is  surrounded  by  an  effective  medium  without  the  presence  of 
the  fluid  envelope.  If  the  porosity  is  relatively  high,  the  prescribed  choice  of  R may  be 
plausible.  If  the  porosity  is  rather  low,  on  the  other  hand,  the  influence  of  neighboring 
elements  may  not  be  represented  by  the  porosity  only.  In  this  case,  the  detailed  structure 
of  the  neighboring  elements  layout  should  have  strong  influence.  Recently,  Dodd  et  al. 

(1995)  recognized  this  aspect  and  adopted  R = ac[l  - S^)]172  / (1-  e)yi  for  two 

dimensional  flow  through  fibrous  media,  where  S(0)  is  the  zero-wave-number  structure 

factor  of  the  cylindrical  element  layout.  The  theoretical  reason  for  this  choice  of  R can  be 

found  in  Dodd  et  al.  (1995).  For  a random  layout  of  cylinders,  the  expression  for  S(0)  has 

been  given  by  Chae,  Ree  and  Ree  (1969)  as  follows: 

()  = [(1  - 1.9682(1  - e)  + 0.9716(1  - s)2  f 

U 1 + 0.0636(1  - e)  - 05446(1  - e)2  - 0.4632(1  - e)3  - 0.1060(1  - s)4  + 0.0087(1  - e)5 

(3-3) 

We  have  found  that  the  three  dimensional  analogue  of  this  choice  of  R is  not  as 
good  as  the  one  determined  by  porosity  matching.  So  the  R calculated  from  the  structure 
factor  will  only  be  applied  to  two-dimensional  flow.  For  the  purpose  of  presenting  the 
results,  we  will  denote  the  cases  of  R=l,  R calculated  by  porosity  matching,  R determined 
from  S(0)  as  EMA-I,  EMA-II  and  EMA-III,  respectively. 
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The  boundary  conditions  for  equations  (3-1)  and  (3-2)  are:  no-slip  condition  at 
r=a,  continuity  of  velocity  and  stress  at  r=aR  and  the  prescribed  velocity  field  u=U  far 
away  from  the  collector  (i.e.,  at  r— »oo). 


The  dimensionless  form  of  the  momentum  equations  are: 


Vp  = V2u 

1 < ? < R 

(3-4) 

Vp  = V2u  - a2u 

? > R 

(3-5) 

where  p = p / [pU  / a ] , 

u = u/U,  a = a/Vk,  and  ? = r/a. 

For  brevity,  the  carets 

indicating  dimensionless  variables  are  dropped  in  the  following. 

The  solutions  to  equations  (3-3)  and  (3-4)  for  the  different  hydrodynamic  cases  of 
flow  through  fibrous  media  are  presented  in  the  following  sections  of  this  Chapter  and  for 
flow  through  granular  media  will  be  given  in  Chapter  4. 


3-2  Classification  of  Random  Packings  of  Cylinders 

According  to  their  hydrodynamic  behaviors,  flow  through  fibrous  media  may  be 
classified  into  the  following  categories:  Flow-I,  All  fiber  (circular  cylinder)  axes  are 
perpendicular  to  the  superficial  velocity.  Fiber  axes  all  lay  in  planes  which  are  parallel  to 
each  other  but  perpendicular  to  the  superficial  fluid  velocity  while  the  fibers  can  have 
random  angles  in  the  planes,  as  shown  in  Figure  3.2(a).  This  case  is  of  great  practical 
importance  since  many  fibrous  media  are  prepared  by  depositing  fibers  onto  a flat  surface 
with  their  axes  almost  completely  parallel  to  the  surface.  A good  example  is  filter  paper 
prepared  by  passing  a slurry  of  fibers  through  a screen  on  which  the  fibers  are  deposited. 
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Another  case  in  which  the  fibers  are  in  a bundle  across  which  the  fluid  flows  is  a special 
case  of  this  category,  as  seen  in  Figure  3-2  (b).  Flow-II  Flow  through  completely 
randomly  packed  fibers.  Fiber  axes  are  completely  randomly  orientated  in  all  possible 
directions.  Figure  3-2  (c).  Flow-Ill,  Planar  packing  of  fibers  with  mean  flow  parallel  to 
the  packing  planes.  Fiber  axes  all  lay  in  planes  parallel  to  each  other  and  to  the  superficial 
fluid  velocity  but  have  completely  random  angles  in  the  planes,  as  seen  in  Figure  3-2(d). 
This  category  is  encountered  with  media  prepared  as  in  the  first  case  of  Flow-I,  but  the 
flow  is  now  in  the  direction  parallel  to  the  planes  on  which  the  fibers  fall.  Flow-IV,  All 
fibers  are  parallel  to  the  mean  flow.  Fibers  are  in  a bundle  and  the  fluid  flows  parallel  to 
their  axes,  as  shown  in  Figure  3-2(e). 


(c)  Flow-II 


U 


(e)  Flow-IV 


Figure  3-2  Hydrodynamic  classification  of  random  packings  of  cylinders 
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3-3  Packing  of  Uniform  Circular  Cylinders  Parallel  to  Mean  Flow 


This  is  corresponding  to  the  case  of  Flow-IV  shown  in  Figure  3-2  (e).  The 
schematic  of  the  EMA  model  for  this  case  is  shown  in  Figure  3-3.  The  dimensionless 
form  of  the  momentum  equations  (3-4)  and  (3-5)  for  the  current  case  are  reduced  to: 


dp  1 d f duz  . D 

dz  r dr  v dr 


(3-6) 


dp  _ 1 d 
dz  r dr 


f du  z 
r — - 
V dr  / 


-a//Uz,  r>R 


(3-7) 


Figure  3-3  Schematic  of  the  EMA  model  for  flow  parallel  to  circular  cylinders 


The  boundary  conditions  include:  (a)  no-slip  condition  at  cylinder  surface  (r=l); 
(b)  fluid  velocity  tends  to  the  superficial  velocity  far  away  from  the  cylinder;  (c)  flow 

field  continues  at  r = R ; (d)  stress  continues  at  r = R . For  constant  -dp/dz(=  af/  ),  the 
solutions  to  equations  (3-6)  and  (3-7)  are  readily  obtained  as: 

2 

u7=- — - + Alnr  + B,  l<r<R 
z 4 dz 


(3-8) 
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uz=CI0(a?,r)  + DK0(<4r)-^A_,  r>R 

dz  a ii 


(3-9) 


where  A,  B,  C,  D are  integral  constants,  I0  and  K0  are  Modified  Bessel  functions  of  the  0th 
order.  The  integral  constants  are  determined  by  using  the  four  boundary  conditions  as 


A = - 


dp 


Raz/K^az/R) 


dz  K0(a//R)  + a/zRKl(a/zR)lnR 


( 1 1 R2  R2  . 

— - — + InR 

\ct//  4 4 2 


B = - 


1 dp 
4 dz 


C = 0 


D = 


dp 


1 1 | R2  R2 

aj,  4 + 4 2 


InR 


^ Rfdp 
j 2 dz 


(3-10) 

(3-11) 

(3-12) 

(3-13) 


dz  K0(azzR)  + azzRK|(azzR)lnR 
In  the  equations  above,  K,  is  the  modified  Bessel  function  of  the  second  kind  of  the  Is 
order. 


The  drag  force  acting  on  the  cylinder,  FD// , can  then  be  calculated  as  follows 


^ZL  = 27tdUz 


dr 


r_,=  2d  — — + IK 
r_1  1 2 dz 


(3-14) 


Simple  force  balance  gives  the  relationship  between  the  drag  and  pressure  gradient  as  the 
following: 


dp  _ 1 - e Fd// 
dl  tc  Z 


(3-15) 


In  the  above  equations,  the  length  coordinate  for  each  cylinder  is  denoted  by  z, 
while  Z is  the  total  length  of  the  cylinders.  The  depth  coordinate  for  a porous  medium 
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bed  is  denoted  by  1,  while  L is  the  depth  of  the  bed.  Substituting  equations  (3-10)  and  (3- 
14)  into  equation  (3-15)  yields  the  following  transcendental  equation  for  determination  of 


a 


//• 


i(i+— 

2 V l-e 


1 R2  R2 


Va// 


h - 

2 4 4 


In  R 


a//RK1(ct//R) 


+ ■ 


R" 


K0(a//R)  + a//RK1(a//R)ln  R 


(3-16) 

The  dimensionless  permeability  K/; , which  is  only  a function  of  the  porosity  e,  is 
then  obtained  once  an  is  solved  from  (3-16): 


K//  = 


//  - 2 
a// 


(3-17) 


Table  3-1  Permeability  equations  by  different  models  for  randomly  packed  uniform 
circular  cylinders  parallel  to  mean  flow 


Model 

(author) 

2 a2 
«//  =T“ 
k// 

Cell  Model 
(Happel,1959) 

a2  8(1_8) 

" 1 — 4e  — 2 ln(l  — e)  — 4(1  — e)2 

Brinkman 
(Spielman  & 
Goren, 
1968) 

1 2K,(a/7) 

1-8  a//K-o(a//) 

EMA  Model 
(present) 

'fl+  1 I 

^ 1 J + R2_R2|nR] 

a//RKl(a//R)  ! 

2I  1-J 

Wn  4 4 2 J 

Kota/yRi  + az/RK^a/yRilnR  2 
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Figure  3-4  Permeability  comparison  between  experimental  data  and  model  predictions 
for  packing  of  uniform  circular  cylinders  parallel  to  mean  flow 


Table  3-1  is  a summary  of  the  permeability  equations  for  Happel’s  model, 
Brinkman’s  model  and  the  EMA  model  obtained  currently. 

Figure  3-3  shows  the  comparison  between  the  permeability  predictions  by  EMA 
and  Happel’s  models  and  the  experimental  data  by  Sullivan  (1942).  Sullivan’s 
measurement  were  carried  out  by  flowing  conditioned  air  parallel  through  bundles  of 
straightened  goat  wool,  Chinese  hair,  glass  wool,  blond  hair  and  copper  wires.  The 
tabulated  measurement  parameters  were  converted  to  the  relations  between  K and  e by 
the  author.  Considering  the  fact  that  it  is  difficult  to  align  the  thin  fibers  perfectly  in 
parallel  to  the  flow  direction  at  low  porosities  and  the  experiments  tend  to  underestimate 
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the  permeability,  the  EMA  models  with  both  choices  of  the  size  of  the  liquid  envelope 
(EMAI  and  EMA  II)  agree  with  the  experimental  data  reasonably  well. 


3-4  Packing  of  Uniform  Circular  Cylinders  Perpendicular  to  Mean  Flow 


The  solutions  to  the  dimensionless  equations  (3-4)  and  (3-5)  for  fluid  flow 
through  a packing  of  uniform  circular  cylinders  perpendicular  to  the  superficial  velocity 
are  readily  obtained  by  the  use  of  stream  functions: 

1 dvj/ 


u = — 


r ae 


Uo 


d\\i 


10  ~ dr 

Using  i|/,  equation  (3-4)  and  (3-5)  are  reduced  to 
V4v|/  = 0 1 < r < R 


(3-18) 

(3-19) 


V>  =a2iV2v1/ 


r < R 


(3-20) 

(3-21) 


The  far  field  condition  ( i.e.,  u-»U  as  r — >oo)  suggests  that 

\|/(r,  0)  = f(r)  sin0  (3-22) 

Using  the  no-slip  condition  for  (3-20)  and  the  far  field  condition  for  (3-21),  the  solution 
to  (3-20)  and  (3-21)  are  given  as 


VR,  = 


M>2  = 


f , n 

f 

r r3^ 

r3  -2r + - C + 

2r  + rlnr-  — - — 

D 

v rJ 

2 2) 

- 

A 

r h Ba  i Ki(a  i r) 

sinO 

r 

sin  9 


(3-23) 


(3-24) 
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Here  K,  is  the  modified  Bessel  function  of  the  first  order.  The  subscript  1 and  2 indicate 
the  fluid  envelope  region  and  the  effective  medium  region,  respectively. 

The  integration  constants  A,  B,  C and  D in  equations  (3-23)  and  (3-24)  are 
determined  using  the  continuity  conditions  at  r=R  (i.e.,  continuity  of  velocity,  normal  and 
shear  stresses).  A set  of  four  linear  algebraic  equations  for  these  constants  are  obtained 
as: 


A - ^[2K,  (Raj  + axRK0(Ra1  )]B  + (l  - R4  )C  + (r4  - R2  = 0 


Ba2RK0(Rax ) + 2(2R  - R3  + R In  R)D  + 4(r3  - r)c  = 2R 


(3-25) 

(3-26) 


( a]  4 ''' 

+ 


R R J 


-B, 


^2a*K0(Ra±)  4ariKl(Ra,1) 


V R 


R-  J 


+ 41  — r + R IC-4I  — + — |D  = -afR 


R 


R 2 


A—r-B 
R3 


4a1K,(Ra1)  2a2K0(RaJA 


a3  K,(RaJ  + ^ 1 + 


R 


(3-27) 


+ 4|  R + |C  - 2RD  = 0 


(3-28) 

Once  the  stream  function  is  determined,  the  pressure  field  is  obtained  from  the 
momentum  equations  (3-4)  and  (3-5)  as 


P.  = 


( 2^ 

8rC- 

4r  + — D 

k r J 

cos0 


(3-29) 


P2  =- 


7 a, 
a\x  + A — 
r 


cos  6 


(3-30) 


The  tangential  and  normal  stresses  acting  on  the  surface  of  cylindrical  element  are 
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CTreir=.  = -4(2C-y)sin0 


c,  , =2(D-2C)cos0 


(3-31) 


(3-32) 


Thus  the  drag  force  per  unit  length,  which  has  been  scaled  by  pU  can  be  obtained  as 

In 

FD1  = |[crIT|r=i  cos#  - cr^,  sin#jd#  = 4;zD  (3-33) 

o L 

From  equations  (3-25)-(3-28),  D is  obtained  as 


D = 


~Q~ 


{[-2«2R2  + 3R4a2  +16R2  -a2]K,(a1R)  + 8R3a1K0(a1R)}  (3-34) 


where 


Q = 


-al  +2a±R2  +alR4  lnR-a^R4  -a^  lnR  + 16R2a±  In  R 
RK,(a±R)  + [4aiR4  In  R +4a2  R2  -3a2  R4  -ai  + 16R2k0(a±R) 


(3-35) 


The  overall  pressure  drop  in  the  direction  of  flow  is  obtained  by  multiplying  the 
drag  force  per  unit  cylinder  length  by  total  length  of  cylinders  per  unit  volume  of  the 
porous  medium,  (1  -s)/7t: 


dp  _ 1 - e 

“ ~ hD± 

dl  7C 


Substituting  equation  (3-33)  into  (3-36)  yields. 


- — = 4D(l-e) 
dl 


(3-36) 


(3-37) 


This  pressure  gradient  should  equal  to  or2  as  obtained  from  equation  (3-5)  with  a uniform 
approaching  (or  superficial)  velocity.  Thus, 


a2=4D(l-s) 


(3-38) 
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This  equation  along  with  equation  (3-34)  represents  an  implicit  equation  for  a± . Once 
aL  is  determined,  the  pressure  drop  (or  gradient)  is  given  by  equation  (3-37). 


Table  3-2  Permeability  expressions  by  different  models  for  randomly  packed  circular 
cylinders  of  uniform  size  perpendicular  to  mean  flow 


Model 

(author) 

2 a! 

Cell  model 
(Happel,1959) 

a2-  40  "e) 

1 1 , „ * (1-e)2 

“ — ~ — ln(l  - s)  H — y 

2 2 l+(l-s)2 

Cell  model 
(Kuwabara,1959) 

2 _ 4(1-8) 

i i i 

ln(l  - e)-  e — (1  — s)2 

4 2 4 

Brinkman’s 
model  (Spielman 
& Goren,1968) 

a2  = 4(1  - e) 

1 2 . ai K,(ai) 

UC  | ~r 

L2  K0(a2)  J 

EMA 

(present) 

fl-L  -4(1  e)  q • 

[-2«fR2  +3RWl  + 16R2  -afjK^a.R) 
+ 8R3aLK0(aLR) 

► 

Similar  procedure  described  above  can  be  applied  to  the  cell  models  and 
Brinkman  model  to  obtain  the  expression  for  permeability.  In  table  3-2,  the  results  for 
these  models  are  given  for  comparison.  As  mentioned  previously,  the  Brinkman’s  model 
(which  was  worked  out  by  Spielman  and  Goren(1968))  is  a special  case  of  the  EMA 


model  with  R=l. 


56 


Figure  3-5  Comparison  of  theoretical  predictions  of  the  dimensionless  permeability  for 
fluid  flow  perpendicular  to  cylinder  axes.  (The  lines  are  predictions  by  different 
theoretical  models,  the  symbols  are  the  results  by  numerical  calculations) 


Typical  results  for  permeability  are  shown  in  Figure  3-4,  where  the  dimensionless 
permeability  K=k/a2  is  given  as  a function  of  porosity.  Also  plotted  were  the  results  of 
numerical  calculations  by  Wang  and  Sangani  (1997),  and  Koch  and  Ladd  (1997),  which 
have  used  different  numerical  techniques  of  recent  development.  Both  of  them  are  based 
numerical  solution  of  the  Stokes  equation  for  the  random  array  of  cylinders.  In  the 
absence  of  experimental  results,  their  numerical  results  may  serve  as  the  basis  for 
comparison  between  various  models.  At  high  porosities,  all  models  agree  well  with  the 
numerical  results.  As  the  porosity  is  decreased,  the  predictions  of  each  model  start  to 
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show  some  differences.  The  EMA-III  shows  the  best  agreement  with  the  numerical 
results  followed  by  EMA-II  and  the  cell  models  in  order.  EMA-I  , on  the  other  hand, 
shows  significant  underestimation  when  the  porosity  is  smaller  than  about  0.7.  When 
porosity  is  smaller  than  about  0.45,  the  EMA-I  does  not  even  have  a solution  for 
permeability.  This  significant  under  prediction  may  be  a consequence  of  not  having  the 
fluid  envelope  which  exists  in  the  other  models.  Considering  the  simplicity  of  these 
models,  the  permeability  predictions  appear  to  be  quite  reasonable  in  comparison  with 
numerical  and  experimental  results. 


3-5  Planar  Packing  of  Uniform  Circular  Cylinders  with  Mean  Flow 
Parallel  to  Packing  Planes 

This  case  is  corresponding  to  Flow-Ill  as  shown  in  Figure  3-2  (d).  If  3 is  the 
angle  between  a cylinder  axis  and  the  direction  of  the  superficial  velocity,  then  the  drag  in 
the  direction  of  flow  is 

Fd  = FD1  sin  & + Fd//  cos  9 (3-39) 

where  FD1  and  FD//  represent  the  drags  contributed  by  the  flow  components  normal  to  and 
parallel  to  the  cylinder  axis,  respectively. 

In  a two  dimensional  planar  random  packing  of  uniform  cylinders,  i.e.,  the  case  in 
which  the  cylinder  axis  all  lay  randomly  in  a series  of  planes  parallel  to  each  other  and 
the  fluid  flows  parallel  to  the  planes,  the  probability  of  finding  a cylinder  lying  at  any 
angle  is  the  same.  The  probability  density  P8(&)  of  finding  a cylinder  lying  at  an  angle  d 
can  be  found  through  the  following  normalization  condition: 
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|p,(J)d*«l  (3-40) 

Here  we  need  only  to  distinguish  the  cylinder  angle  in  0-90°  range,  because 
cylinders  with  angles  of  0,  S+7t/2,  d+n  and  0+3ti/2  have  the  same  effect  to  flow  in  a 
completely  random  packing. 


Figure  3-6  Drags  acting  on  an  oblique  cylinder 


In  a two  dimensional  planar  random  packing,  P9  does  not  depend  on  3. 
Therefore,  we  have 


pe  =~  (3-41) 

n 

fd  =Jo/2[FDi(d)sin&  + FD//(9)cosS]p0d&  (3-42) 


From  equation  (3-36),  we  have 
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7ta± 


71 


FD±(3)  = -j sinS  = 

1-8  (1-8)K1 


sin  9 


(3-43) 


where  sind  term  reflects  the  fact  that  the  only  the  perpendicular  component  of  the 
superficial  velocity  should  be  used  to  calculate  the  perpendicular  drag,  which  is  scaled  by 
the  pU.  Fd//  is  calculated  from  equation  (3-15): 


7UX// 


71 


fd//(S)  = - cos»  = 

1-8  (l-e)K/; 


cos  & 


(3-44) 


Substituting  equations  (3-43)  and  (3-44)  into  equation  (3-42)  and  carrying  out  the 
integration  gives: 


Fd  = 


71 


2(1-8) 


1 1 
+ • 


Because 


dp  _ 1 - 8 

dT“_ 


71 


D 


where  1 is  in  the  direction  of  the  superficial  velocity.  Therefore 


K ~ 2 


1 1 


vK±  K;/, 


or 


a2  =^(al  + <*//) 


(3-45) 


(3-46) 


(3-47) 


(3-48) 


The  a values  turns  out  to  be  in  a very  simple  relations  with  a //  and  a± , mathematical 


average  of  the  two. 
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3-6  Completely  Randomly  Packed  Uniform  Circular  Cylinders 

In  the  case  of  a completely  random  packing,  i.e.,  the  case  in  which  cylinder  axes 
are  orientated  in  all  possible  directions,  the  probability  of  finding  the  number  of  cylinders 
which  lay  in  a certain  solid  angle  (see  Appendix  A,  “Introduction  to  Solid  Angle”)  is 
proportional  to  the  steradian  of  the  solid  angle.  If  we  denote  P(l)  as  the  probability  density 
of  finding  a cylinder  lying  in  solid  angle  oo,  then  the  normalization  condition  gives, 

Jo2Xdco  = l (3-49) 

For  an  isotropic  three  dimensional  packing,  Pm  does  not  change  with  to,  therefore  one 
obtains 

Pco  (3-50) 

271 

If  the  flow  direction  is  in  the  negative  z direction  as  shown  in  Figure  A,  the  angle 
<|>  where  the  cylinder  axis  lies  does  not  have  any  effect  on  the  drag  acting  on  the  cylinder, 
so  we  need  only  to  find  the  probability  density  of  finding  a cylinder  whose  cone  angle  is 
0,  P8(S).  The  normalization  condition  gives, 

JJi20P8(9)dS=l  (3-51) 

Considering  dco  = sin  9dSd(j> , comparing  equations  (3-49)  and  (3-51)  yields 

P3(S)  = sin0  (3-52) 

Equation  (3-52)  shows  that  the  spatial  cylinder  axis  distribution  in  terms  of  the  cone 
angle  0 is  not  uniform  in  a 3-D  random  packing. 
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The  drag  force  in  the  direction  of  the  fluid  flow  is  then  obtained  similarly  by 
averaging  equation  (3-39)  with  the  weighting  function  P9(9): 

Fd  =Jo/2[FD1(»)sin8  + FD//(»)cos8]P90)d9  (3-53) 


Figure  3-7  Comparison  of  permeability  for  different  hydrodynamic  cases  of  flow  though 
fibrous  media 


From  equations  (3-36)  and  (3-15),  we  have: 


fd±(9)=  1 sin  9 

1 - e 

(3-54) 

Fd//(»)  = 7^»s8 
1 -e 

(3-55) 

Substituting  equations  (3-54)  and  (3-55)  into  equation  (3-53),  we  get 
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Fd  = \q2  sin2  9 + a2  cos2  sjsin  9d9 


(3-56) 


Therefore,  the  a and  K are  obtained  for  the  current  flow  problem  as 


(3-57) 


1 2 1 


(3-58) 


+ 


K 3K±  3K„ 


Figure  3-6  shows  the  comparison  of  permeabilities  for  the  four  hydrodynamic 


cases.  At  the  same  porosity,  the  case  with  flow  normal  to  the  cylinders  results  the  highest 
pressure  drop  while  the  case  with  flow  parallel  to  all  the  cylinders  gives  the  lowest.  The 
permeability  predictions  by  EMA  for  various  flow  cases  through  fibrous  media  are  shown 
in  Figure  3-7  along  with  the  available  experimental  data.  Jackson  and  James  (1986) 


media.  The  earliest  work  was  that  by  Carman  (1938).  In  the  following  year,  Wiggins  et 
al.  (1939)  determined  the  permeability  of  a wide  variety  of  fibrous  media:  glass  rods, 
glass  fiber,  copper  wire  and  Celanese  yarn.  Each  medium  was  a random  packing  of  fibers 
in  a tube,  and  several  test  liquids  were  used.  Sullivan  (1942)  measured  permeability  of 
various  wools  and  hairs.  The  most  important  detail  of  his  work  is  that  the  fibers  were 
packed  in  a tube  such  that  they  were  somehow  aligned  with  the  flow.  Brown’s 
measurement  (1950)  was  for  flow  of  dry  air  through  banks  of  heat  exchanger  tubes 
oriented  normal  to  the  flow.  Davis  (1952)  presented  a lot  of  experimental  data  for  the 


provided  a very  good  collection  of  experimental  data  for  permeability  through  fibrous 


permeability  of  flow  through  packings  of  glass  wool,  glass  rods,  cotton  wool,  camel  hair. 
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kapok,  merino  wool,  rayon,  down,  and  a mixture  of  glass  wool  and  copper  wire.  A few 
years  later,  Chen  (1955)  gathered  data  from  various  experiments  on  filter  mats,  in  which 
the  fibers  were  randomly  oriented  in  planes  normal  to  the  flow.  Some  further  data  on 
filter  mats  were  provided  by  Wheat  (1963)  and  Labrecque  (1968).  Ingmanson  et  al. 
(1959)  measured  the  flow  of  air  through  glass,  nylon  and  paper  fibers  which  were 
compressed  as  the  pressure  drop  across  the  medium  increased.  White  (1960)  obtained  the 
permeability  of  distilled  water  through  various  concentrations  of  cross  linked  acrylamide 
polymer  gels  which  can  be  considered  as  fibrous  media  composed  of  polymer  chains. 
Kirsch  and  Fuchs  (1966)  tested  different  arrangement  of  regular  arrays  of  parallel  fibers 
normal  to  the  flow.  Collagen  is  a protein  which  forms  fibrous  structures  found 
throughout  the  human  body.  The  earliest  measurements  of  collagen  permeability  with 
sufficient  information  to  determine  the  porosity  were  reported  by  Stenzel  et  al.  (1971), 
older  data  are  available  but  the  porosity  information  was  not  available.  Visvanadham  et 
al.  (1978)  reported  the  similar  tests.  Kostornov  and  Shevchuck  (1977)  measured 
permeability  of  water  and  alcohol  through  discs  of  compressed  20%Cr-80%Ni  alloy 
fibers,  50  pm  in  diameter  and  3 mm  long.  Jackson  and  James  (1982)  studied  the 
permeability  of  entangled  hyaluronic  acid  polymer  chains  in  solution.  The  most  recent 
experiments  were  reported  by  Rahli  et  al.  (1996),  who  measured  permeability  of  random 
packing  of  metal  wires  of  different  aspect  ratios  (length/diameter).  Their  permeability 
data  were  lower  than  most  other  experimental  results,  which  was  believed  to  be  caused  by 
the  end  effect  of  the  short  fibers.  Generally,  the  comparison  made  in  Figure  3-7  shows 
that  the  EMA  model  agrees  with  most  of  the  experimental  data  very  well. 
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Figure  3-8  is  a comparison  of  permeability  predictions  of  EMA  and  the  widely 
used  Happel’s  model.  In  view  of  the  comparisons  made  in  Figure  3-7  with  the  big 
quantity  of  experimental  data,  the  EMA  model  proves  to  be  a better  one  than  the 
Happel’s,  which  under-estimates  the  differences  in  flow  patterns. 


3-7  Randomly  Packed  Non-uniform  Circular  Cylinders 


In  the  case  of  a porous  medium  composed  of  poly-disperse  elements,  we  can 
apply  the  analysis  as  we  did  in  Sections  3-3  to  3-6  to  each  of  the  elements.  If  we  scale 
the  length  parameters  by  the  size  of  the  corresponding  element  each  time,  we  can  easily 
see  that  we  would  get  the  same  solutions  for  all  the  elements,  i.e.,  the  same  six  integral 
constants  which  are  only  functions  of  the  porosity  and  element  size.  The  dimensional 
drag  for  each  of  the  elements  can  then  be  obtained  and  the  pressure  drop  calculated. 

From  Section  3-4,  we  know  that  the  drag  acting  on  a cylinder  of  radius  a sitting  in 
a fibrous  medium  of  porosity  e in  the  case  of  mean  flow  perpendicular  to  the  cylinder 
axes  is  given  by 

Fd  = 4rtpD(R,  -£=)UZ  (3-59) 

Vk 


when  a fluid  flows  through  the  medium  with  a superficial  velocity  of  U , where  Z is  the 


length  of  the  cylinder  and  D(R,— =)  is  given  by  equation  (3-34) . 

Vk 
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Figure  3-8  Comparison  of  EMA  predictions  of  permeability  of  fibrous  media  with 
available  experimental  data 


For  a packing  of  polydisperse  cylinders  with  a size-length  distribution  density  of 
Pa(a,Z),  the  number  of  cylinders  with  size  in  the  range  between  a and  a+da  and  length 
between  Z and  Z+dZ  is  then  given  by 

dn(a)  = NPa  (a,  Z)dadZ  (3-60) 

where  N is  the  total  number  of  cylinders  in  the  medium  of  volume  V.  The  total  volume 


of  these  cylinders  is 

Vs  = f0N  Tta2  Zdn(a)  = «N  JJ°  JJ*  Pa  (a,  Z)Za2dadZ 


(3-61) 
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Figure  3-9  Comparison  of  predictions  of  permeability  by  EMA  and  Happel’s  models 


Because  there  are  N cylinders  contained  in  the  porous  medium  of  volume  V , 


i Vs 

1 - E = — = 71  — 

V IvJ 


Jo°  Jo°  Za2Pa(a,Z)dadZ 


(3-62) 


or 

(3-63) 

If  the  lengths  of  the  cylinders  are  much  bigger  than  their  diameters,  the  drag  experienced 
on  the  ends  of  each  of  the  cylinders  can  then  be  neglected.  The  total  drag  acting  on  all 
the  cylinders  in  volume  V is  then  obtained  from  equation  (3-59), 
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Ft  =4tihUNJ®J®  ZD(R, -j=)Pa  (a,  Z)dadZ  (3-64) 

v k 

The  drag  force  is  balanced  by  the  pressure  drop  in  a column  of  porous  medium  of  length 
L and  cross  sectional  area  of  A 

~TE=Ia=4’vUV,0”,0”°  ZD<R'-/=>p.<a’z>dadz  (3-«) 


Substituting  equation  (3-63)  into  equation  (3.65),  we  get  the  following  expression  for  the 
pressure  drop  in  a porous  medium  composed  of  polydisperse  cylinders  perpendicular  to 
the  mean  flow: 


~~  = 4pU(l  - s) 


Ho”  ZD 


R, 


yfk- 


Pa(a,Z)dadZ 


Jo°  Jo°  Za2Pa(a,Z)dadZ 


(3-66) 


For  the  special  case  of  uniform  cylinders,  equation  (3-66)  is  reduced  to  equation 
(3-37)  discussed  in  Section  3-4.  Therefore,  the  equation  for  k is  then 


f = 4(1-8)- 
k 


JoT  ZD 


R, 


Vk. 


Pa(a,Z)dadZ 


lo°  C Za  2 Pa  ( a,  Z )dadZ 


(3-67) 


This  equation  is  the  generalized  permeability  equation  which  can  be  used  to  calculate 
permeability  of  a polydisperse  packing  of  cylinders  perpendicular  to  the  superficial 
velocity.  In  a polydisperse  packing,  the  thickness  of  the  clear  liquid  region  around  each 
solid  element  may  be  considered  the  same.  The  thickness  5 can  be  estimated  through  the 


overall  packing  porosity. 
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J J a3Pa(a,Z)dadZ 


1-8  = 


0 0 


f J (a  + 6)3  Pa(a, Z)dadZ 

o o 


(3-68) 


R in  the  permeability  equation  is  then  obtained  by  R=  (a+8)/a.  If  the  size  distribution  is 
not  very  broad,  equation  (3-67)  may  be  approximated  by 


f*4(l-e)- 
k 


D 


R(am)- 


Jo*  Jo”  ZPa(a,Z)dadZ 


Jo°  Jo°  Za2  Pa  (a,  Z)dadZ 


(3-69) 


where  the  cylinder  size  a in  function  D[R,a  / Vk]has  been  replaced  by  a middle  value  a„,. 
Equation  (3-69))  can  also  be  written  as 


^-*4(1-8)d[  R(am),  am 


Vk- 


(3-70) 


where 


ae  = 


lo°  Jo°  Za2Pa(a,Z)dadZ^  '2 


J„T  ZPa(a,Z)dadZ 


(3-71) 


is  an  average  size,  which  may  be  a good  estimate  for  am,  i.e., 


m e 


(3-72) 


Therefore,  equation  (3-69)  is  then  reduced  to  the  following 


w 4(1  - e)of  R(ae  ),—£= 
k V Vk 


(3-73) 


which  has  the  same  form  as  equation  (3-38).  This  implies  that  equation  (3-71)  is  the 
definition  of  an  “equivalent  radius”  of  a polydisperse  packing,  the  physical  meaning  of 
which  can  be  understood  as  the  radius  of  the  uniform  cylinders  packed  to  give  the  same 
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permeability  as  the  medium  composed  of  polydisperse  cylinders  of  the  same  porosity.  It 
is  also  noted  that  R(ae)=l/(l-g)l/3,  it  is  only  a function  of  the  porosity.  It  can  also  be 
easily  shown  that  extending  the  results  for  other  types  of  packing  described  in  Sections  3- 
3 to  3-6  to  the  corresponding  polydisperse  cases  gives  similar  conclusion.  The 
significance  of  equation  (3-71)  is  that  all  the  equations  derived  before  for  the  cases  of 
uniform  cylinders  are  also  valid  for  polydisperse  cylinders  if  we  replace  the  radius  of  the 
cylinders  with  ae. 


3-8  Summary 

The  Effective  medium  approximation,  which  assumes  a model  system  in  which  a 
packing  element  (a  single  fiber  in  the  fibrous  medium  and  a single  sphere  in  the  granular 
medium)  is  surrounded  by  a fluid  envelope  and  an  effective-medium  beyond  the 
envelope,  is  applied  to  flow  through  fibrous  media  to  predict  the  permeability  of  four 
different  flow  cases:  Flow-I,  the  mean  velocity  is  perpendicular  to  the  axes  of  all  fibers; 
Flow-II,  the  mean  flow  is  parallel  to  the  packing  planes  in  which  all  fibers  are  randomly 
oriented;  Flow-Ill,  flow  through  completely  randomly  packed  fibers  and  Flow-IV  in 
which  the  superficial  velocity  is  parallel  to  the  axes  of  all  the  fibers.  The  new  model 
integrated  the  important  features  of  both  the  cell  models  and  Brinkman’s  model.  The 
analytical  expressions  for  the  permeability  of  flow  through  fibrous  media  are  obtained. 
Comparisons  with  the  experiment  and  the  strictly  numerical  calculated  results  showed 
that  the  EMA  based  on  the  porosity  matching  and  structure  factor  are  better  flow  models 
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than  the  mostly  used  cell  and  Brinkman’s  models.  The  EM  A analysis  is  also  extended  to 
the  fibrous  media  with  polydisperse  elements. 

The  criteria  for  the  derived  equations  to  be  valid  may  include  the  following.  The 
aspect  ratio  (length/diameter)  of  the  fibers  should  be  sufficiently  large  so  that  the  effect  of 
the  fiber  ends  could  be  neglected.  The  Reynolds  number  of  the  flow  should  be  small 
enough  in  order  to  have  the  Stokes  equation  to  hold.  It  is  generally  accepted  that  a 
Reynolds  number  of  1 is  a critical  one  below  which  the  flow  is  considered  to  be  creeping 
flow.  But  engineering  practice  shows  that  this  critical  number  could  be  extended  to  10. 
Another  criterion  is  that  the  Knudsen  number  (length  of  the  mean  free  path  / fiber 
diameter)  should  be  less  than  0.01  if  the  fluid  is  gas,  because  the  no-slip  condition  has 
been  used  at  the  fiber  surfaces. 


CHAPTER  4 

EFFECTIVE  MEDIUM  APPROXIMATION  FOR  GRANULAR  MEDIA 


4-1  Packing  of  Uniform  Spheres 


As  pointed  out  previously,  the  important  feature  of  the  EMA  model  is  that  it 
reflects  the  physical  situation  for  the  flow  around  the  representative  element  more 
realistically  than  the  cell  models  which  solve  only  equation  (3-1).  The  flow  field  around 
each  of  the  solid  elements  which  constitute  the  packed  medium  collectively  is  affected  by 
the  neighboring  elements  which  are  represented  by  an  effective  medium  in  the  current 
model.  While  each  solid  element  is  surrounded  by  an  effective  medium  which  has  the 
same  permeability  as  the  entire  porous  medium,  it  contributes  to  that  property  at  the  same 
time.  Thus,  the  permeability  of  the  porous  medium  and  the  local  flow  field  around  each 
solid  element  are  coupled  together.  The  EMA  model  accounts  for  this  coupling,  whereas 
the  cell  models  do  not.  This  more  realistic  feature  of  the  EMA  model,  however,  results  in 
a somewhat  more  complicated  solution  procedure  although  an  analytic  solution  is  still 
obtained. 

The  simple  geometry  of  the  model  system  allows  analytic  solution  to  these 
equations  by  the  use  of  stream  function: 


u = - 


1 di|/  2 f (r) 


r2  sinG  50 


cos0 


(4-1) 
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1 cty 
rsin0  dr 


1 df(r) 
r dr 


sinG 


where 


(4-2) 


v|/  = f(r)  sin20  (4-3) 

and  u has  been  scaled  to  dimensionless  by  U,  the  superficial  velocity  of  the  fluid  in  the 
packed  bed.  All  the  length  parameters  have  been  scaled  by  the  radius  of  the  spheres,  a. 

Far  away  from  the  particle,  the  velocity  field  approaches  the  uniform  superficial 
velocity.  Thus,  in  the  effective-medium  region  (r  > R),  the  radial  dependence  of  the 
stream  function  f(r)  is  given  as 


1 


f,(r)  = — A + a + - 


a2r 


e^B- 


r/ 


(4-4) 


Here  a2  is  a2/k  which  is  the  inverse  of  dimensionless  permeability.  In  the  fluid  envelope 
region  (1  < r < R),  on  the  other  hand,  the  no-slip  condition  at  the  particle  surface  gives  the 
function  f(r)  as 


f2(r)  = 


f 4 2 i 

r r 1 
C 20+U~30rJ 


C + 


f _2  , i X 

r r 1 

<T_2  + 6r, 


D 


(4-5) 


The  subscripts  1 and  2 in  Equations  (4-4)  and  (4-5)  represent  the  effective-medium  region 
and  the  fluid  envelope  region,  respectively. 

The  four  constants  A,  B,  C,  and  D appearing  in  Equations  (4-4)  and  (4-5)  are 
determined  by  the  matching  conditions  at  the  surface  of  the  fluid  envelope  (i.e.,  at  r = R), 
which  represent  the  continuity  of  velocity  and  stresses.  These  matching  conditions 
provide  a set  of  four  linear  equations  which  gives  an  analytic  expression  for  A,  B,  C,  and 


D in  terms  of  a and  R.  The  four  linear  equations  are: 
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(4-7) 


(4-8) 


(4-9) 


Using  the  velocity  described  above  (i.e.,  Equations  (4-1)  and  (4-2)  along  with  (4-4) 
and  (4-5)),  the  pressure  field  ( scaled  by  pU/a  ) in  each  region  is  determined  from  the 
momentum  equations  (3-1)  and  (3-2)  as  follows: 


. 2 1 
P.  =1  -°rr  + ^7A, 


COS0 


P 2 =lrC  + ^TDy 


cos0 


(4-10) 


(4-11) 


Since  the  velocity  and  the  pressure  fields  are  now  known,  the  drag  force  ( scaled  by 
pUa)  on  individual  particle  can  be  determined  as 
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Fd  = |27i|ar0|r  a sin0 - arr[  acos0jsin0d0 


(4-12) 


where  the  normal  and  shear  stresses  ( scaled  by  pU/a)  at  the  particle  surface  are  given  as 


a_  =(C  + D)cos0 

ri  I r=a 


’re]r= 


( 1 

— C - D sin0 
V2 


(4-13) 

(4-14) 


Thus,  the  drag  force  on  a single  spherical  particle  is 


Fd  = -4nD 


(4-15) 


Assuming  that  the  particles  are  uniformly  distributed  in  the  packed  bed,  the  pressure  drop 
over  a unit  distance  in  the  flow  direction  is  a multiple  of  FD  and  the  number  density  of  the 
particle  3(1-s)/4tt.  Therefore,  the  pressure  gradient  in  the  flow  direction  1 is  given  as 


dp  3(1  - e) 
dl  4ti 


Fd=3(1-8)D 


(4-16) 


From  the  Darcy’s  law  (or  Equation  (3-2)),  the  pressure  gradient  for  a uniform  flow  is 
1/K.  Therefore,  the  expression  for  K is  obtained  as 


a2=^  = -3(1-e)D 

The  solution  of  the  Equations  (4-6)-(4-9)  gives  the  following  expression  for  D. 


(4-17) 


D = -6 


6R6a3  + 21a2R5  - 5R4a3  +45aR4 
- 5R3a2  + 45R 3 -a3R  -a2 


(4-18) 


where 


Q = 1 80R3  +24RV  -45R4a2  -9a2  +4a3  + 1 80aR4  +4R6a3 
+ 10(aR)3  -9Ra3  -9a3R5  +30(aR)2  - 1 80R3a 


(4-19) 


75 


Equation  (15)  and  (16)  uniquely  determine  the  value  of  the  permeability  K once  the 
porosity  s (hence  R)  is  specified. 


4-2  Packing  of  Non-uniform  Spheres 


In  this  section,  we  consider  the  case  of  a porous  medium  composed  of 
polydisperse  spheres  with  a size  distribution  density  of  P(a).  For  the  convenience  of 
presentation,  all  the  variables  in  this  section  are  dimensional.  The  number  of  particles 
with  their  radii  in  the  range  of  a and  a + da  is  then  given  as 

dn(a)  = NP(a)da  (4-20) 

Here  N is  the  total  number  of  particles  in  the  porous  medium  of  volume  V . The  total 
volume  of  these  N particles  is 

Vs=f0"|’ta3‘ln(a)  = yNM3  (4-21) 

where 


M3  = J”a3P(a)da  (4-22) 

is  the  third  moment  of  the  particle  size  distribution.  Since  Vs/V  represents  the  volume 
fraction  occupied  by  the  solid  particles,  the  porosity  of  the  packed  bed  is  then  given  as 


1 -e 


Zs 

V 


4ti(NV, 

M3 

3 \VJ  i 


(4-23) 


Therefore  the  particle  number  density  of  the  porous  medium  is 

N _ 3(1  -s) 

V ~~  4tiM3 


(4-24) 
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Since  the  drag  force  on  an  individual  particle  of  radius  a sitting  in  a porous 
medium  of  porosity  s is  given  as  Equation  (13),  the  total  drag  force  on  all  the  particles  in 
a unit  volume  is  obtained  as 


Ft  = -47tpU^Jjj°  aD(R,  a )P(a)da 


(4-25) 


Since  this  total  drag  is  equivalent  to  the  pressure  gradient  in  the  porous  medium, 
the  following  expression  for  the  permeability  of  a porous  medium  composed  of 
polydisperse  spheres  is  obtained  as 
1 3(1-8), 


Mi 


-J“aD(R,a)P(a)da 


(4-26) 


The  thickness  of  the  clear  liquid  envelope,  8,  is  estimated  form  the  overall  porosity  of  the 
packing. 


1 -8  = 


Ja’Pa(a)da 
_o 

}(a  + S)3Pa(a)da 


(4-27) 


R in  the  permeability  equation  (4-26)  is  then  obtained  by  R=1.028(a+8)/a.  The  factor 
1.028  will  be  discussed  in  the  next  section.  For  the  special  case  of  a bidisperse  packing, 
equation  (4-27)  is  reduced  to: 


( 8 'l 

3 

f,  51 

1 + — 

xi 

+ 

1 + — 

V a J 

( a 2J 

1 


(4-28) 


where  x(  is  the  volume  fraction  of  particle  size  aj. 


For  the  special  case  of  uniform  spheres  of  size  a,  P(a)  = 8(a),  here  8(a)  is  the  delta 
function.  Thus  M3  = a3  from  Equation  (4-22),  and  Equation  (4-26)  is  reduced  to  Equation 
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(4-17)  described  in  the  last  section.  Thus,  the  permeability  k of  a packed  bed  filled  with 
polydisperse  spherical  particles  is  represented  by  an  implicit  equation  (4-26)  along  with 
Equation  (4-18).  Although  this  equation  is  complicated  algebraically,  the  permeability  k 
is  rather  easily  determined  by  a simple  numerical  scheme  once  the  particle  size 
distribution  P(a)  and  the  porosity  e are  specified. 


Porosity  e 


Figure  4-1  Comparisons  of  permeability  predictions  for  packing  of  uniform  spheres 
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4-3  Comparison  of  Permeability  Predictions 


In  this  section,  the  permeability  predicted  by  the  present  model  is  compared  with 
those  obtained  through  pure  numerical  simulations  and  predicted  by  other  models 
including  the  Kozeny-Carman  correlation  for  a random  packing  of  uniform  spheres  and 
non-uniform  spheres.  The  numerical  results  for  packing  of  uniform  spheres  by  Ladd 
(1990),  Kang  and  Sangani  (1994)  and  Kim  and  Russel  (1985)  are  almost  identical. 
Essentially  all  these  investigators  calculated  the  permeability  from  a multipole-moment 
expansion  of  force  density  on  the  surface  of  solid  particles.  The  porosity  range  reported 
by  Phillips  et  al.  (1988)  and  Kim  and  Russel  (1985)  is  >0.5,  Ladd  (1990)  is  >0.55,  and 
Kang  and  Sangani  (1994)  is  > 0.39.  The  Kozeny-Carman  correlation  is  represented  by 
the  following  equation: 


(4-29) 


This  equation  was  initially  developed  by  Kozeny  based  on  the  capillary  tube  model,  and 
later,  Carman  adjusted  the  numerical  constant  empirically  from  16  to  45  to  provide  better 
agreement  with  experimental  results  (Carman,  1956;  Tien,  1989).  Thus,  Equation  (4-29) 
is  a semi-empirical  correlation.  This  Kozeny-Carman  equation  has  been  found  to  provide 
good  agreement  with  numerous  experimental  observations  for  the  porosity  range  of  0.26 
to  0.8  which  is  achieved  in  packed  beds  and  sedimentation  experiments  (Happel  and 


Brenner,  1965).  Carman  (1956)  also  found  that  Equation  (4-29)  could  be  applied  to 
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mixtures  of  various  particle  sizes  if  the  hydraulic  radius  was  used  in  place  of  the  particle 
radius. 


Figure  4-2  Comparison  of  permeability  prediction  by  EMA  model  and  experimental  data 
for  a bidisperse  particle  system.  (The  radii  of  the  big  and  small  particles  are  147  nm  and 
1 03  nm,  respectively) 


For  the  simplicity  of  presentation,  the  EMA  models  with  choices  of  R=l,  R=  1/(1- 
s)l/3  and  R=[(l  -S(0))/(l-e)]13  are  denoted  as  EMA  I,  EMA  II  and  EMA  III,  respectively. 
S(0)  is  the  zero  wavenumber  structure  factor  of  the  sphere  packing.  The  zero 
wavenumber  structure  factor  S(0)  is  defined  by 


S(0)  = J[P(r|0)  - n]dr 


(4-30) 


80 


where  P(r|0)  is  the  probability  for  finding  a particle  with  its  center  in  the  vicinity  of  r 
given  a test  particle  at  origin,  and  n is  the  number  density  of  particles.  S(0)  for  randomly 
distributed  spheres  is  given  by  the  well-known  Carnahan-Starling  formula 


S(0)  = 


l + 4(l-e)  + 4(l-s)2  -4(l-e)3  +(1-e)4 


(4-31) 


0 0.2  0.4  0.6  0.8  1 

Volume  faction  of  small  particles 


Figure  4-3  Comparison  of  permeability  prediction  by  EMA  model  and  experimental  data 
for  a bidisperse  particle  system  (The  radii  of  the  big  and  small  particles  are  605  nm  and 
230  nm,  respectively). 


Figure  4-1  indicates  that  the  cell  model  of  Happel  shows  a very  good  agreement 
with  the  numerical  simulation  results.  Considering  its  simplicity,  the  Happel’s  cell  model 
appears  to  be  an  excellent  model  for  the  prediction  of  the  permeability  of  packed  beds 
filled  with  uniform  spheres.  We  may  also  note  that  the  Brinkman’s  model  is  poorer  than 
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either  model.  When  the  porosity  is  smaller  than  about  0.4,  this  model  cannot  even 
provide  a permeability  prediction  as  the  solution  becomes  singular.  EMA  II  and  EMA  III 
agree  with  the  numerical  results  very  well  for  porosities  higher  than  about  0.5,  while 
EMA  III  is  a little  better  than  EMA  II. 


Figure  4-4  Logarithmic  normal  distributions  of  packing  spheres 

The  permeability  of  a packed  bed  is  directly  associated  the  drag  force  on  the 
packing  particles  (or  the  energy  dissipation  at  the  surface  of  the  particles).  When  the 
porosity  is  high,  the  physical  contacts  between  packing  particles  may  be  sparse,  and  the 
simplifying  assumption  that  a spherical  fluid  envelope  encapsulates  the  particle  may  not 
introduce  a serious  error  since  the  drag  force  is  mostly  determined  by  the  flow  field  in  the 
immediate  vicinity  of  the  particles.  With  decreasing  porosity,  however,  the  number  of 
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physical  contacts  between  packing  particles  increases.  Consequently,  the  assumption 
which  redistributes  the  fluid  uniformly  around  the  particles  may  become  more  unrealistic 
and  is  likely  to  introduce  a more  significant  error.  This  argument  may  explain  the  small 
discrepancy  observed  with  the  EMA  model  at  a low  porosity  in  comparison  with  the 
Kozeny-Carman  correlation  and  numerical  results. 


Porosity,  £ 

Figure  4-5  Permeability  predictions  by  EMA  and  Kozeny-Carman  correlation  for 
polydisperse  packing  (Lines:  EMA  predictions  (Eqn.  4-26);  Symbols:  Kozeny-Carman 
correlation) 


Although  it  is  a pure  empiricism,  a minor  adjustment  of  increasing  the  radius  of 
liquid  envelope  R by  a mere  2.8  percent  from  1/(1  -s)l/3  to  1.028/(l-s)1/3  is  found  to 
improve  the  agreement  with  the  numerical  results  providing  an  excellent  match  for  the 
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whole  porosity  range.  The  Kozeny-Carman  correlation  clearly  overestimates  the 
permeability  when  the  porosity  is  higher  than  about  0.6. 

The  major  contribution  of  the  present  work  may  be  the  extension  of  the  EM  A 
model  to  the  case  of  polydisperse  packings.  As  far  as  we  know,  the  present  approach 
described  in  Section  3 is  the  first  attempt  which  fully  incorporates  the  particle  size 
distribution  for  the  permeability  prediction.  The  permeability  of  a packed  bed  with 
polydisperse  particles  has  been  usually  estimated  using  the  Kozeny-Carman  equation 
with  the  reciprocal  mean  particle  size.  The  reciprocal  mean  particle  size  is  the  average 
size  which  corresponds  to  the  same  total  specific  surface  area  of  packing  particles.  This 
approach  has  been  known  to  provide  a good  prediction  for  the  permeability  of  spherical 
as  well  as  non-spherical  packings  if  their  particle  size  distribution  is  not  very  broad 
(Brenner,  1961;  Standish  and  McGregor,  1978;  Perry,  1984). 

In  Figure  4-2  and  Figure  4-3,  comparisons  are  made  between  the  permeability 
predictions  of  the  EMA  model  using  the  adjusted  R (i.e.,  R=l. 028/(1  -s)l/3),  we  will  call  it 
EMA  IV,  and  the  experimental  data  by  Thies-Weesie  and  Philipse  (1994),  who 
investigated  flow  of  cyclohexane  through  bidisperse  silica  spheres  for  various  mixing 
ratios.  The  radii  of  the  big  and  small  particles  are  103nm  and  147nm  in  Figure  4-2,  and 
605  nm  and  230  nm  in  Figure  4-3,  respectively.  In  the  calculation  the  measured  porosity 
data  as  measured  by  the  same  authors  were  used.  The  measured  porosity  data  have  error 
bars,  so  we  have  two  predicted  curves  for  each  case.  Also  shown  in  the  same  figures  are 
the  prediction  of  Kozeny-Carman  correlation.  It  is  shown  that  both  the  EMA  IV  model 
and  the  Kozeny-Carman  correlation  give  very  good  predictions. 
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In  Figure  4-5,  a comparison  is  made  between  the  predictions  of  the  EMA  model 
using  the  adjusted  R and  the  Kozeny-Carman  equation  for  the  particle  size  distributions 
given  in  Figure  4-4.  Although  the  logarithmic  normal  distribution  has  been  used  in  the 
present  study  because  of  its  wide  use  in  particle  sizing,  any  type  of  particle  size 
distribution  can  be  used.  The  size  distributions  shown  in  Figure  4 have  the  average  size 
(number  averaged)  of  1 with  an  increasing  variance  from  0.00,  corresponding  to  uniform 
spheres,  to  0.284  which  represents  a very  broad  distribution.  The  agreement  between  the 
two  predictions  is  reasonably  good  at  low  porosities,  and  the  Kozeny-Carman  correlation 
overestimates  the  permeability  at  higher  porosities  as  is  the  case  in  monodisperse 
packings. 

In  the  sense  of  predicting  the  permeability,  the  EMA  model  is  more  complicated 
than  the  Kozeny-Carman  correlation.  The  advantage  of  EMA  model  over  the  Kozeny- 
Carman  correlation  is  that  it  gives  a flow  field  around  each  of  the  solid  elements,  which  is 
very  important  for  many  heat  and  mass  transfer  processes  in  packed  beds. 

4-4  Summary 

The  permeability  of  packed  beds  filled  with  polydisperse  spheres  has  been  studied 
using  the  Effective  Medium  Approximation  (EMA).  The  EMA  assumes  a model  system 
in  which  a packing  element  (or  particle)  is  surrounded  by  a fluid  envelope  and  an 
effective-medium  beyond  the  envelope.  Since  the  local  flow  field  around  each  packing 
element  which  is  influenced  by  the  neighboring  elements  contributes  collectively  to  the 
permeability  of  the  packed  bed,  the  local  flow  and  the  permeability  are  coupled  together. 
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The  EMA  model  reflects  this  physical  situation  of  coupling  more  realistically  than  other 
existing  models.  Furthermore,  unlike  other  models,  this  model  is  capable  of  fully 
incorporating  the  packing  size  distribution  in  predicting  the  permeability  of  packed  beds. 
Although  the  coupling  between  the  local  flow  field  around  an  individual  packing  element 
and  the  permeability  of  the  collective  medium  (i.e.,  packed  bed)  makes  the  solution 
procedure  complicated,  the  simple  geometry  of  the  model  unit  enables  us  to  obtain  an 
analytic  expression  for  the  permeability  of  the  packed  bed  as  a function  of  the  packing 
size  distribution  and  the  porosity.  The  prediction  of  the  present  model  shows  a good 
agreement  with  the  numerical  results  in  the  whole  range  of  porosity  and  the  Kozeny- 
Carman  correlation  at  low  porosities. 

In  modeling  the  effective  heat  or  mass  transfer  coefficient  in  porous  media  (or 
packed  beds),  a representative  flow  field  for  a model  system  is  needed.  While  the 
Kozeny-Carman  correlation  provides  an  excellent  prediction  for  the  permeability  of  a 
porous  medium  at  low  porosities,  a representative  flow  field  is  not  given  for  these 
purposes.  The  cell  models  and  the  EMA  model,  on  the  other  hand,  provide  a 
representative  flow  field  for  the  prediction  of  the  effective  transfer  coefficients.  Despite 
the  simplicity  of  these  models,  the  permeability  prediction  by  them  appears  to  be  good 
probably  due  to  the  fact  that  the  flow  field  in  the  immediate  vicinity  of  the  packing 
element  is  most  important  in  determining  the  drag  force  on  the  packing  elements  (hence 
the  permeability). 

When  the  convection  is  dominant  over  diffusion,  only  the  flow  field  in  the 
immediate  vicinity  of  the  packing  elements  is  important.  Thus,  the  cell  models  or  the 
EMA  model  are  expected  to  provide  a reasonable  prediction  for  the  effective  transfer 


86 


coefficients.  Application  of  the  EMA  model  for  that  purpose  is  our  major  objective,  and 
the  present  chapter  is  a part  of  that  effort.  Since  the  EMA  model  portrays  the  flow 
situation  in  a packed  bed  more  realistically  than  the  cell  models,  we  expect  that  the 
prediction  of  effective  heat  or  mass  transfer  coefficient  by  the  EMA  model  will  be  better 
that  those  by  other  models,  and  the  results  will  be  presented  in  Chapter  6. 


CHAPTER  5 

PARTICLE  DEPOSITION  IN  FIBROUS  AND  GRANULAR  MEDIA 
5-1  Particle  Deposition  on  Circular  Cylinders  Parallel  to  Mean  Flow 

In  Sections  1-3  and  1-4,  it  has  been  shown  that  under  widely  applicable 
conditions,  particle  deposition  on  solid  surfaces  in  the  presence  of  van  der  Waals,  electric 
double  layer,  and  hydrodynamic  interactions  can  be  treated  by  assuming  convective  mass 
transfer  in  the  fluid  and  a first  order  reaction  at  the  collector  surface.  In  order  to  solve  the 
convective  diffusion  equation,  a representative  flow  field  should  be  available.  In  the 
current  case,  i.e.,  liquid  flows  parallel  to  a cylinder  surface,  the  mass  transfer  boundary 
layer  develops  in  the  flow  direction.  In  a laminar  flow  outside  a circular  cylinder,  the 
mass  transfer  boundary  layer  in  the  down  stream  can  be  quite  thick  and  is  limited  by  the 
presence  of  neighboring  cylinders.  We  have  shown  in  the  foregoing  chapters  that  only  a 
flow  field  close  to  the  cylinder  surface  is  predictable  by  using  the  effective  medium 
approximation.  Therefore  we  can  only  solve  the  convective  mass  transfer  equation  in  the 
end  region,  where  mass  transfer  boundary  layer  has  not  fully  developed.  Fortunately,  this 
limitation  in  analysis  does  not  really  pose  a serious  problem,  because  the  deposition  rate 
in  the  down  stream  is  usually  a few  magnitude  lower  than  the  end  region  as  shown  later 
in  this  section.  That  is  to  say,  particle  deposition  occurs  mostly  at  the  very  end  region,  to 
which  the  remaining  part  of  this  section  will  be  dedicated. 
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If  the  Peclet  number  aU/£>  has  a value  much  higher  than  unity,  the  diffusion 
boundary  layer  will  be  sufficiently  thin  in  the  end  region,  so  that  the  curvature  effect  and 
the  tangential  diffusion  are  negligible.  Under  these  conditions,  the  convective  diffusion 
equation  assumes  the  following  form: 


3c 

u7  — = 


32c 


dz  dr 2 


The  dimensionless  form  of  equation  (5-1)  is 


(5-1) 


u , 


3c  1 32c 


3z  Pe  3f2 


(5-2) 


where  uz  = uz  / U,r  = r / a,c  = c / c0,  c0  is  the  bulk  particle  concentration.  For  brevity, 
the  carets  indicating  dimensionless  variables  are  dropped  in  the  following. 

The  above  problem  is  the  so-called  extended  Graetz  problem,  for  which  an 
infinite  series  solution  can  be  obtained  (Bowen  et  al.  1976).  When  the  Peclet  number  is 
large,  a simpler  approximate  solution  which  is  asymptotically  correct  for  the  entrance 
region  can  be  found.  In  the  entrance  region,  the  particle  concentration  boundary  layer  is 
still  developing.  The  simplest  approximate  form  is  obtained  when  this  boundary  layer 
lies  entirely  within  the  region  near  the  collector  surface  where  the  velocity  profile  can  be 
replaced  by  its  tangent  line  at  the  surface.  Thus  the  problem  is  reduced  to  the  one  of 
particle  deposition  on  a surface  in  an  infinite  medium  undergoing  a simple  shear  flow. 
Equation  (5-2)  is  then  reduced  to 

3c  1 32c 

Twy ~ T~~  2 

3z  Pe  dy 


(5-3) 
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where  y=r-l . Let  t wy  = 2Y , and  X = (2  / 9z)l/3  Y , equation  (5-3)  is  transformed  to  the 
following  ordinary  differential  equation 


d2c 


+ 3X2  — = 0 


dXz  dX 

The  corresponding  boundary  conditions  are 
c(X->oo)=l 


(5-4) 


c(X  = 0)  - 


f Z>oo  \ 

1/3 

' dc  A 

vaK0y 

v9zJ 

,dX7 

= 0 


(5-5) 

(5-6) 


x=o 


The  solution  of  equation  (5-4)  satisfying  the  boundary  conditions  (5-5)  and  (5-6)  is 


V 2 A 1/3 


c(X)  = 


v9z  7 


a,  ' 

vaK0y 


-) 

9z/ 


1/3 


(5-7) 


+ 


r(t) 


or  c(z,Y)  = 


8 

f-1 

\9z  ) 

1/3 


+ Jq2/9z)  Y e-C 


v aKjjj , 


f2 

<9z 


1/3 


(5-8) 


r(j) 


The  particle  deposition  flux  is  then  readily  obtained  from  equation  (5-8) 

'2aV/V‘ 


Sh  = 


Ja 

VC  A 


9zJ 


ax, 


UpU, 


( V ^ 


vaK^y 


1/3 


-)  + 

V9z/ 


r(t) 


(5-9) 


All  the  variables  except  Sh  in  equation  (5-9)  are  dimensional.  The  shear  stress  term  can 
be  calculated  from  equation  (3-8)  to  be 
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z^a 


Figure  5-1 . Mass  transfer  coefficient  as  a function  of  the  distance  from  the  ends  of  a 
bundle  of  parallel  cylinders  when  fluid  flows  parallel  to  them 


axw_a2  Ra3K|(aR)  1 1 + R2  R2 

pU  2 K0(aR)  + aRK|(aR)lnR^a2  4 4 2 

(5-10) 

Shown  in  figure  5-1  is  the  Sh  as  a function  of  dimensionless  distance  from  the  end 
of  the  cylinders.  The  packing  with  a lower  porosity  gives  a little  higher  particle 
deposition  rate.  Particle  deposition  is  significant  only  in  the  very  entrance  region  where 
the  diffusion  boundary  layer  is  still  developing.  Even  at  a position  z/a=0.001  and  Da=oo 
(the  surface  reaction  is  instantaneous),  the  Sh  is  extremely  small  compared  to  the  case  of 
flow  perpendicular  to  the  cylinders  where  the  Sh  is  usually  in  the  order  of  102  for  most 
particle  deposition  processes  (see  Section  5-2).  This  treatment  does  not  consider  the 
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particle  deposition  on  the  end  surfaces  of  the  fibers,  which  may  be  neglected  if  the  fiber 
aspect  ratio  (length/diameter)  is  big. 


Table  5-1  The  flow  field  parameter  Ac  for  various  models 


Model  (author) 


Cell  model 
(Happel,  1959) 


A =- 
c 2 


2 2 l + (l-E) 


Cell  model 
(Kuwabara,  1959) 


Ac  ~ 2^-4  ~ 2 ln(1  ~ e)  - 8 ~ 4(1  ~ e)2 


Brinkman’s  model 
(Spielman  & Goren,1968) 


aK  (a)  a2  2 , , Ma) 

A _ = — 1 , where  = 2a  + 4a  ■ 


2K0(a) 


1 - s 


K0(a) 


Ac  = 


aR 

2Q 


-4a2R2(l  + ln  R)  + 3R4a2  +16R2  +a2 


+ 8aR 


R 


-l]K0(aR) 


K,(aR) 


EMA  model 
(present) 


where 


Q = [-a3+2aR2+a3R4lnR-a3R4-a3lnR  + 16R2alnR] 
RK,(aR)  + [4a2R4  In  R + 4a2R2  -3a2R4  - a2  +16R2jK0(aR) 


5-2  Particle  Deposition  on  Circular  Cylinders  Perpendicular  to  Mean  Flow 

As  discussed  in  Section  1-3,  a linearized  form  of  the  flow  field  u is  sufficient  for 
the  calculation  of  concentration  distribution  since  the  convective  diffusion  equation  is 
applied  within  the  diffusion  boundary  layer.  For  the  case  in  which  the  flow  direction  is 
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perpendicular  to  the  axes  of  the  cylindrical  collectors,  the  linearized  velocity  field  is 
represented  by  the  following  stream  function: 

y = -2Ac(r  - 1)2  sin0  (5-11) 

Here  Ac  is  the  flow  field  parameter  which  is  unique  to  each  of  the  porous  medium 
models.  For  the  EMA  model,  Ac  can  be  obtained  from  the  analysis  shown  in  Section  3-4 
to  be: 


A„  = 


aR  | 
2QI 


-4a2R2(l  + lnR)  + 3R\U  +16R^  +a 


4 2 


K,(aR)  + 8aR 


R 


-l]K0(aR)J 


(5-12) 


Where  Q is  given  by  equation  (3-35).  The  flow  field  parameter  Ac  for  other  models  is 
summarized  for  comparison  in  Table  5-1. 

Using  the  IFBL  approximation,  Spielman  and  Friedlander  (1974)  obtained  the 
dimensionless  particle  deposition  rate  as: 


Sh  = 0.73 1 A |./3  Pe  !/3y  e (p  c ) (1-38) 

In  Figure  5-2,  the  mass  transfer  coefficients  predicted  by  the  EMA  and  the  cell 
models  are  presented  as  a function  of  the  packing  porosity.  Also  plotted  here  are  the 
recent  numerical  results  of  Wang  and  Sangani  (1997).  Wang  and  Sangani’s  work  is  not 
for  a particle  deposition  process  but  for  a heat  transfer  problem  in  which  the  overall  heat 
transfer  coefficient  has  been  calculated  numerically  for  a flow  through  an  array  of 
randomly  positioned  cylinders  which  are  held  at  constant  temperature.  Since  Wang  and 
Sangani’s  results  do  not  involve  any  assumptions  other  than  a high  Peclet  number,  they 
may  serve  as  a basis  to  evaluate  the  effectiveness  of  various  models.  Figure  5-2  indicates 
that  the  EMA  model  provides  very  good  predictions  for  the  mass  transfer  coefficient  for  a 
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broad  range  of  porosity.  The  cell  models,  on  the  other  hand,  appear  to  overestimate  the 
Sh  by  about  20%. 


Figure  5-2.  Sh  predicted  by  various  models  for  packing  of  fibers  perpendicular  to 
superficial  velocity 


5-3  Particle  Deposition  on  Planar  and  Completely  Randomly  Packed  Circular  Cylinders 

The  particle  deposition  in  the  cases  of  Flow-Ill  and  Flow-IV  is  treated  in  this 
section  based  on  the  results  from  Sections  5-1  and  5-2.  The  fluid  flow  relative  to  each  of 
the  fibers  can  always  be  considered  to  be  the  vector  summation  of  a perpendicular 
component  and  a parallel  component  relative  to  the  fiber  axis  section  wise  due  to  the 
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linearity  of  the  Stokes  equation  and  the  convective  diffusion  equation.  As  discussed  in 
Section  5-1  and  5-2,  the  contribution  of  the  parallel  component  to  the  total  particle 
deposition  is  significant  only  in  the  very  end  region,  where  the  diffusion  boundary  layer 
still  undergoes  development.  The  fiber  aspect  ratio  is  rather  high  in  most  applications. 
The  particle  deposition  on  the  ends  of  the  fibers  and  the  particle  deposition  due  to  the 
parallel  component  of  the  fluid  flow  may  be  neglected. 

If  we  denote  the  angle  between  the  axis  of  a section  of  fiber  and  the  superficial 
velocity  by  $,  then  the  perpendicular  component  of  the  superficial  velocity  contribution  is 

U j_  = U sin  & (5-13) 

The  Sh  for  this  section  of  fiber  can  then  be  calculated  by  ( from  equation  1-38) 

Sh(9)  = 0.73  lAj./3Pel/3yc(Pc)sin 1/3  8 (5-14) 

The  overall  Sh  is  then  obtained  as  a probability  average 

Sh  = 0.731Aj./3Pe1/3yc(Pe)j,Q/2P9({))sin1/3  8d3  (5-14) 

where  P9($)is  the  probability  density  of  finding  a fiber  laying  at  an  angle  relative  to  the 
superficial  velocity,  which  has  been  obtained  for  Flow-Ill  and  Flow-IV  in  Sections  3-5 
and  3-6,  respectively.  Equation  (5-14)  may  also  be  written  in  the  following  form: 

Sh  = 0.731A"V,3yc(Pc)x  (5-15) 

where  / = Jq  Pg (9  )sin 11  9d9,  which  can  be  considered  as  a correction  factor  for  the 
packing  structure.  For  a Two-dimensional  planar  random  packing, 

Z-fjj'W*  9dS  = = 0.8235  (5-16) 
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and  in  the  case  of  completely  randomly  packed  fibers, 


X = Jo/2sin4/3  »d&  = 


r(I)r(!) 

2r(f) 


= 0.9107 


(5-17) 


Table  5-2  The  flow  field  parameter  As  for  various  models 

Model  (author) 

As 

HappeTs  cell  model 
(Spielman  & FitzPatik,  1973) 

2f  1 — (1  — e)5/3l 

A J 

' s 2 - 3(1  - e),/3  + 3(1  - e)5/3  - 2(1  - e)2 

Brinkman’s  model 
(Present) 

As  = 1 + a 

EMA  model 
(Present) 

2 [~-30R4a3  -30R3a2  +90R3  +90aR4  +42a2R5 

As  = — 

Q[-15a2R2  +12R6a3  +15a3R3  +3a2  +3a3R 
where 

Q = 1 80R3  + 24R5a2  - 45R4a2  - 9a2  + 4a3  + 1 80aR4  + 4R6a 
+ 10(aR)3  -9Ra3  -9a3R5  +30(aR)2  -180R3a 

5-4  Particle  Deposition  on  Randomly  Packed  Non-uniform  Circular  Cylinders 


The  particle  deposition  on  randomly  packed  non-uniform  cylinders  can  be  treated 
readily  based  on  the  results  of  Sections  5-2  and  5-3  if  the  size  distribution  Pa(a)  and  the 
surface  area  fraction  available  for  each  size,  PA(a),  or  the  length  fraction  for  each  size  are 
known.  The  surface  averaged  flux  can  then  be  calculated  by 


J = 0.731Xc0U  1/3Z>2/3J0”a;/3 


v 


£.-?=■)  a 2/3yc(s,a)Pa(a)PA(a)da 


VkJ 


(5-18) 
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here,  the  packing  structure  correction  factor  x has  been  assumed  to  be  independent  of  the 
cylinder  size.  The  integral  in  equation  (5-18)  can  only  be  calculated  through  numerical 
procedures  due  to  the  complicated  functions  At(a)  and  yc(a).  This  case  will  not  be 
discussed  further  because  there  are  no  experimental  data  available  yet  in  the  literature  to 
compare  with. 


Porosity 


Figure  5-3.  Sh  predicted  by  various  models  for  granular  media 


5-5  Particle  Deposition  on  Randomly  Packed  Uniform  Spheres 


A representative  linearized  velocity  field  for  flow  passing  through  granular  media 


is  represented  by  the  following  stream  function: 
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H/  = -|As(r-l)2sin2e  (5-19) 

Here  As  is  the  flow  field  parameter  which  is  unique  to  each  of  the  porous  medium 
models.  For  the  EMA  model,  As  can  be  obtained  from  the  analysis  shown  in  Section  4-1 
to  be: 


As  = 


-30R4a3  -30RV  +90R3  +90aR4  +42a2R5 
- 1 5a  2 R 2 + 12R6a3+15a3R3  +3a2  +3a3R 


(5-20) 


where  Q is  given  by  equation  (4-19).  The  flow  field  parameter  As  for  other  models  is 
summarized  for  comparison  in  Table  5-2. 

Using  the  IFBL  approximation,  Spielman  and  Friedlander  (1974)  obtained  the 
dimensionless  particle  deposition  rate  as: 


Sh  = 0.653A'/3Pel/3ys(Ps) 


(1-38) 


In  Figure  5-3,  the  mass  transfer  coefficients  predicted  by  the  EMA  and  the  cell 
models  are  presented  as  a function  of  the  packing  porosity.  Also  plotted  here  are  the 
recent  numerical  results  of  Sorensen  and  Stewart  (1994a,  1974b).  Sorensen  and  Stewart’s 
determined  the  heat  and  mass  transfer  coefficients  for  simple  cubic  and  face-centered 
cubic  arrays  at  their  closest  packing  porosities  of  s=0.4764  and  0.2595,  respectively.  The 
modified  EMA-II  model,  i.e.,  the  EMA  model  where  the  size  of  the  liquid  envelope  has 
been  modified  to  accommodate  the  matching  of  the  permeability  predictions  with  the 
experimental  results,  R=1.028/(l-s)l/3,  gives  exact  prediction  for  the  face-centered  array, 
although  it  overestimates  the  mass  transfer  rate  for  simple  cubic  array  by  about  12%.  As 
mentioned  before,  the  EMA  model  is  essentially  developed  for  random  arrays.  It  can 
neither  distinguish  the  exact  packing  structure  nor  the  flow  direction  relative  to  the 
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structure.  Also  shown  in  the  same  plot  is  the  experimental  results  by  Wilson  and 
Geankoplis  (1966).  The  experimental  data  appeared  to  be  lower  than  the  EMA  prediction 
by  about  12%.  It  is  still  too  early  to  conclude  which  model  gives  better  prediction  at  this 
stage.  This  requires  well  conditioned  experiments  or  numerical  algorithms  for  randomly 
packed  spheres. 


5-6  Particle  Deposition  on  Randomly  Packed  Non-uniform  Spheres 


The  particle  deposition  on  randomly  packed  non-uniform  spheres  can  be  analyzed 
based  on  the  results  of  Section  5-5  if  the  size  distribution  Pa(a)  is  known.  The  surface 
averaged  flux  is 


J = 0.635c0Ui/3£>2/3J^  aJ./3 


( a ' 
e’~7r 

v Vk. 


a-2/3yc(8,a)Pa(a)PA(a)da 


(5-21) 


where  PA(a)  is  the  fraction  of  surface  area  available  for  collector  size  a,  which  can  be 
obtained  from  the  size  distribution  density,  Pa(a). 


PA(a)  = 


a2Pa(a) 

J“a2Pa(a)da 


(5-22) 


5-7  Summary 

Effective  medium  approximation  is  applied  to  predict  particle  deposition  rate  in 
fibrous  and  granular  media.  The  analytical  solutions  for  the  particle  deposition  are 
obtained.  Comparison  with  available  experimental  data  shows  that  the  EMA  model  for 
fibrous  media  is  superior  to  the  other  existing  models.  Well  conditioned  experiments  are 
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needed  for  checking  the  validity  of  the  EMA  model  for  granular  media.  The  EMA  model 
is  also  extended  to  treat  particle  deposition  in  polydisperse  fibrous  and  granular  media. 

The  treatment  here  is  not  limited  to  the  deposition  of  colloidal  particles.  It  can  be 
applied  to  other  heat  and  mass  transfer  processes  in  laminar  flow  through  fibrous  and 
granular  media. 


CHAPTER  6 

INFLUENCES  OF  VARIOUS  SYSTEM  PROPERTIES  AND  OPERATING 
CONDITIONS  ON  PARTICLE  DEPOSITION  RATE 


6-1  Influence  of  System  Properties 

Representative  model  predictions  for  various  physical,  chemical  properties  and 
the  operating  condition  on  the  rate  of  deposition  of  colloidal  particles  will  be  presented  in 
this  section.  This  will  be  carried  out  for  deposition  onto  cylindrical  collectors 
perpendicular  to  the  superficial  velocity  based  on  the  effective  medium  approximation.  It 
should  be  emphasized  that  the  general  dependence  of  the  particle  deposition  rate  on 
various  parameters  for  other  media  is  similar  to  the  one  discussed  here. 

Figure  6-1  shows  the  effects  of  Da  and  Pe  on  Sh.  For  the  case  where  the 
Damkohler  number  (Da=Koa/D00)  » 1,  the  particle  deposition  is  convective  diffusion 
limited  and  the  filter  coefficient  is  independent  of  K^.  This  is  equivalent  to  assuming  that 
a perfect  particle  sink  exists  on  the  collector  surface.  When  Da«l,  the  particle 
deposition  is  then  surface  interaction  limited.  In  most  cases  the  particle  deposition  is 
almost  impossible  in  this  instance.  The  Figure  6-1  can  be  considered  as  a master  plot 
which  correlates  the  mass  transfer  coefficient  with  all  the  affecting  parameters.  The 
effects  of  these  individual  parameters  will  be  discussed  below. 
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Figure  6- 1 Influence  of  Da  and  Pe  on  Sh 
6-1-1  Surface  Potentials 

Figure  6-2  illustrates  the  effect  of  particle  and  collector  surface  potentials  on  the 
mass  transfer  coefficient  at  various  superficial  velocities.  The  surface  potentials  affect 
the  magnitude  of  the  double  layer  interaction.  From  the  expressions  of  the  double  layer 
interaction,  it  is  noticed  that  it  is  the  product  of  the  surface  potentials  of  the  collector  and 
particle  that  comes  into  effect.  As  observed,  the  mass  transfer  coefficient  at  a given 
superficial  velocity  decreases  dramatically  as  the  product  of  the  surface  potentials 
increases  above  a certain  value.  Also  plotted  in  Figure  6-2  is  the  corresponding  height 
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of  the  repulsive  energy  barrier.  The  energy  barrier  is  very  sensitive  to  the  surface 
potentials,  and  it  markedly  increases  as  the  product  of  the  surface  potentials  increases 
above  a critical  value.  The  height  of  the  energy  barrier  corresponding  to  this  critical 
value  is  about  5 kBT.  As  the  height  of  the  energy  barrier  increases,  the  mass  transfer 
coefficient  decreases  rapidly  and  beyond  an  energy  barrier  of  about  20  kBT,  the  particle 
deposition  becomes  virtually  impossible.  It  is  also  found  from  our  numerous 
simulations  that  the  particle  deposition  behavior  interpreted  in  terms  of  height  of  energy 
barrier  is  almost  independent  of  the  system,  i.e.,  the  critical  energy  barriers  of  5 kBT  and 
20  kBT  are  not  dependent  on  the  system  parameters  such  as  particle  size,  Hamaker 
constant  and  electrolyte  concentration.  A critical  height  of  repulsive  energy  barrier  of 
about  20  kBT  has  long  been  observed  by  many  studies  on  particle  aggregation 
(Elemelech  et  al.  1995).  When  the  repulsive  energy  barrier  is  higher  than  20  kBT  a 
colloidal  suspension  is  usually  very  stable.  In  fact  particle  deposition  and  aggregation 
are  essentially  the  same  physical  phenomenon  governed  by  the  same  natural  law. 
Particle  deposition  can  be  considered  to  be  a special  case  of  particle  aggregation  in 
which  small  colloidal  particles  aggregate  with  relatively  huge  particles — collectors. 

6-1-2  Ionic  Concentration 

Ionic  concentration  affects  the  mass  transfer  coefficient  through  the  double  layer 
interaction.  The  effect  is  shown  in  Figure  6-3.  It  is  observed  that  the  mass  transfer 
coefficient  increases  as  the  electrolyte  concentration  increases.  This  is  caused  by  the 
decrease  of  the  range  and  magnitude  of  the  repulsive  double  layer  interaction,  which  is 
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usually  called  the  compression  of  the  electrical  diffusion  double  layer.  The  results  shown 
in  Figure  6-3  also  demonstrate  that  the  effect  of  electrolyte  concentration  is  much  more 
profound  for  bigger  particles  than  for  smaller  particles. 
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Figure  6-2  Influence  of  particle  and  collector  surface  potentials  on  mass  transfer 
coefficient  (Particle  diameter=0.5  pm,  electrolyte  concentration=0.01M,  Hamaker 
constant=8.5xl0’21  J,  1-1  electrolyte,  temperature  = 298  K) 


Shown  in  Figure  6-4  is  the  effect  of  particle  size  on  the  relative  magnitudes  of  the 
van  der  Waals  and  double  layer  interactions.  For  particles  much  bigger  than  the  van  der 
Waals  interaction  range  (e.g.,  10  nm),  both  the  double  layer  and  van  der  Waals  interaction 
energies  are  approximately  proportional  to  the  size  of  the  particle.  If  the  particle  size  is  in 
the  van  der  Waals  interaction  range,  the  van  der  Waals  interaction  energy  becomes 
relatively  independent  of  the  particle  size  while  the  double  layer  interaction  energy  still 
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decreases  as  the  particle  size  is  reduced.  Therefore  the  van  der  Waals  force  usually 
dominates  over  the  double  layer  force  for  extremely  small  particles.  In  this  instance  the 
mass  transfer  coefficient  becomes  relatively  insensitive  to  the  changes  in  surface 
potentials  and  electrolyte  concentration. 
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Figure  6-3  Influence  of  electrolyte  concentration  on  mass  transfer  coefficient 
(Superficial  velocity=lxl0'3  m/s;  surface  potentials  (pc=cp(=0.02V,  Hamaker  constant 
=8.5x1 0'21  J;  1-1  electrolyte;  temperature  = 298  K) 

6-1-3  Particle  Size 

The  results  shown  in  Figure  6-3  also  demonstrate  the  effect  of  particle  size  on  the 
mass  transfer  coefficient.  Particle  size  comes  into  the  picture  through  the  particle 
diffusion  coefficient  and  the  size  effect  on  interaction  energies.  Decreasing  the  particle 
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size  results  in  a higher  particle  diffusivity  and  a lower  repulsive  energy  barrier,  and 
therefore  a higher  particle  deposition  rate. 

The  current  research  focuses  only  on  deposition  of  colloidal  particles,  in  which 
the  diffusion  and  convection  are  the  main  mechanisms  for  particles  to  be  transported 
toward  the  collector  surfaces.  For  deposition  of  bigger  particles,  other  mechanisms  such 
as  particle  inertia,  interception  may  play  major  roles.  It  has  been  recognized  that  particles 
of  sizes  around  1pm  are  the  most  difficult  to  be  captured.  For  particles  bigger  than  1pm, 
the  deposition  rate  may  increase  as  the  their  size  is  increased. 


Separation  Distance  (m) 


Figure  6-4  Influence  of  particle  size  on  double  layer  and  van  der  Waals  interactions 
(Electrolyte  concentration  = 0.00 1M,  1-1  electrolyte;  Hamaker  constants. 5x  10-21  J; 
surface  potentials  (pc=cpp=0.02V;  temperature  = 298  K) 
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6-1-4  Hamaker  Constant 

The  van  der  Waals  interaction  energy  is  proportional  to  the  Hamaker  constant  of 
the  system.  As  the  Hamaker  constant  is  increased,  the  mass  transfer  coefficient  increases 
due  to  the  stronger  attractive  force.  This  is  shown  in  Figure  6-5.  Once  the  van  der  Waals 
force  becomes  strong  enough,  corresponding  to  an  energy  barrier  of  about  5 kBT,  the 
particle  deposition  is  limited  by  the  convection  and  diffusion  of  particles  in  the  fluid 
medium  and  the  mass  transfer  coefficient  becomes  independent  of  the  Hamaker  constant. 


Hamaker  Constant  (J) 


Figure  6-5  Influence  of  Hamaker  constant  on  mass  transfer  coefficient  (Particle  diameter 
=0.5pm;  electrolyte  concentration  = 0.01  M,  1-1  electrolyte;  surface  potentials  cpc=  cpp= 
0.02V;  temperature  = 298  K) 
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Superficial  Velocity,  U (m/s) 


Figure  6-6  Influence  of  superficial  velocity  on  mass  transfer  coefficient 


6-2  Influence  of  Operating  Conditions 

The  operating  conditions  may  include  the  superficial  velocity  and  suspension  pH 
because  these  two  parameters  are  often  adjusted  in  the  process  to  achieve  optimum 
operation.  The  suspension  pH  affects  the  particle  deposition  rate  by  changing  the  surface 
potentials  of  the  particle  and  collector  and  the  electrolyte  concentration.  The  correlation 
of  surface  potentials  as  a function  of  pH  are  usually  obtained  through  experiments  for 
different  materials  and  tabulated  in  literature.  The  ionic  concentration  can  be  calculated 


readily  from  the  value  of  pH  and  concentrations  of  other  ions.  Therefore  the  pH  comes 
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into  the  picture  through  physical  properties  of  the  system,  which  has  been  discussed  in 
Section  6-1.  The  effect  of  superficial  velocity  is  discussed  below. 

The  superficial  velocity  is  a measure  of  the  degree  of  convection.  Figure  6-6 
illustrates  the  effect  of  superficial  velocity  on  mass  transfer  coefficient.  When  the 
particle  deposition  is  convective  diffusion  limited  (the  height  of  repulsive  energy  barrier 
is  lower  than  about  5 kBT),  the  mass  transfer  coefficient  is  proportional  to  1/3  power  of 
the  superficial  velocity.  As  the  role  of  surface  repulsive  force  becomes  more  important, 
the  contribution  of  convection  decreases  gradually.  Although  the  particle  deposition  rate 
per  unit  area  of  the  filter  surface  increases  with  superficial  velocity,  the  filter  coefficient, 
which  is  a measure  of  the  speed  of  concentration  decay  decreases  with  increasing 
superficial  velocity  because  the  residence  time  of  the  particles  in  the  filter  is  reduced. 

6-3  Effect  of  Hydrodynamic  Interaction 

As  a particle  approaches  very  close  to  the  surface,  it  becomes  increasingly 
difficult  for  the  liquid  between  them  to  drain  out  of  the  gap  and  this  tends  to  contribute  a 
repulsive  force  to  the  total  interactions  between  the  particle  and  the  surface.  This  effect  is 
usually  termed  as  the  hydrodynamic  interaction. 

One  approach  to  incorporate  the  contribution  of  the  hydrodynamic  interaction  into 
the  particle  deposition  mechanism  is  to  lump  the  effect  to  an  effective  diffusion 
coefficient,  Ve  (Ruckenstein  and  Prieve  1 974).  The  particle  mobility,  «,  defined  as  the 
velocity  produced  by  applying  a unit  force  on  it,  and  the  effective  diffusion  coefficient, 
Ve,  are  related  by  the  Nernst-Einstein  equation: 
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Z>,  ="kT  (6-1) 

The  mobility  of  Brownian  particles  suspended  in  liquids  can  be  determined  by  the  well 
known  Stokes  equation: 


1 

671  pa 


(6-2) 


where  p is  the  viscosity  of  the  liquid  and  a is  the  particle  radius.  The  Stokes  equation 
provides  good  agreement  with  experimental  results  obtained  at  small  Reynolds  numbers, 
i.e.  Re  <1.  Some  of  the  theoretical  treatments  of  this  subject  has  grown  out  of  the  early 
work  of  Stokes.  An  extension  to  account  for  the  inertial  effect  has  been  made  by  Oseen 
(1910). 
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(6-3) 


in  witch  U is  the  velocity  of  the  particle  relative  to  the  liquid,  p the  density  of  the  liquid, 
which  has  been  shown  to  be  valid  for  Re=2aUp/p  <5  (Schlichting,  1960). 

When  the  particle  is  very  near  to  the  solid  surface,  the  mobility  is  decreased  as  a 
result  of  friction  between  the  surface  and  fluid  which  increases  the  force  required  to  push 
the  liquid  out  of  the  path  of  the  approaching  particle.  By  convention,  the  Stokes  equation 
is  modified  by  a dimensionless  factor  g(h),  so  that: 


1 

m 6g(h)7ipa  ^ ^ 

The  hydrodynamic  factor,  g(h),  is  a function  of  the  separation  distance,  h . Some 
analytical  expressions  corresponding  to  the  translational  motion  of  single  particles  normal 
to  the  surface  are  give  in  Table  6-1.  The  rotational  motion  of  particles  is  often  considered 
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to  have  negligible  effects  on  particle  transport  leading  to  deposition.  In  the  range  of 
separation  distances  when  hydrodynamic  effect  needs  to  be  considered  (h«a),  the 
motion  of  Brownian  particles  normal  to  the  surface  is  dominant  because  of  the  strong  van 
der  Waals  and  Double  layer  interactions,  and  the  negligible  fluid  velocity  relative  to  the 
surface.  So  the  hydrodynamic  factors  corresponding  to  the  motion  of  particles  parallel  to 
the  surface  are  not  listed  in  Table  6-1. 


Table  6-1  Hydrodynamic  correction  factors 
Correction  factor  g(h)  Validity  Source 
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Fig.  6-7  shows  the  hydrodynamic  correction  factor  as  a function  of  the  separation 
distance.  Since  1 > 1/  g(h)  > 0,  the  effective  diffusion  coefficient  is  always  smaller  in 
the  region  near  the  surface  than  in  the  bulk  and  it  vanishes  when  the  particle  touches  the 


surface. 


Ill 


After  considering  the  hydrodynamic  effect,  the  equation  for  determining  the 
virtual  chemical  reaction  constant  K now  becomes: 

K«=2„{f[g(h)e«kT-l]dh}‘'  (1-35) 

The  influence  of  hydrodynamic  effect  on  the  virtual  chemical  reaction  constant  is  plotted 
in  Figure6-8.  It  can  be  seen  that  the  hydrodynamic  effect  can  decrease  the  reaction 
constant  by  a factor  of  100.  Figure  6-9  shows  the  influence  of  the  hydrodynamic  effect 
mass  transfer  coefficient.  The  hydrodynamic  effect  is  limited  in  a very  narrow  region 
close  to  the  surface,  it  can  be  considered  to  be  a kind  of  surface  “force”.  So,  it  does  not 
have  effect  on  particle  deposition  if  the  process  is  convective  diffusion  limited. 


h/a 


Figure6-7  Hydrodynamic  retardation  factor  as  a function  of  separation  distance 
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Product  of  Surface  Potentials,  cpccpp  (V2) 

Figure  6-8  Influence  of  hydrodynamic  effect  on  virtual  reaction  constant 


Product  of  Surface  potentials,tpcipp(V2) 


Figure  6-9  Influence  of  hydrodynamic  effect  on  mass  transfer  coefficient  (Particle 
diameter  = 0.25pm;  granular  medium  with  collector  diameter  of  1mm;  porosity  = 0.6; 
Hamaker  constant  = 8.5xl0'21  J,  temperature  = 298K;  electrolyte  concentration  = 0.01M, 
1-1  electrolyte) 
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6-4  Summary 

Effects  of  various  system  properties  and  operating  conditions  on  deposition  of 
colloidal  particles  have  been  discussed.  The  physical  or  chemical  conditions  include  the 
properties  which  affect  the  magnitude  of  double  layer  interaction,  the  electrolyte 
concentration  and  surface  potentials,  and  the  property  which  affects  the  van  der  Waals 
interaction,  the  Hamaker  constant.  Increasing  the  ionic  concentration  and  decreasing  the 
product  of  the  surface  potentials  can  enhance  the  particle  deposition  rate  by  reducing  the 
magnitude  of  the  repulsive  double  layer  interaction.  The  attractive  van  der  Waals 
interaction  is  proportional  to  the  Hamaker  constant,  which  depends  on  the  properties  of 
the  particle  and  collector  and  also  the  suspending  medium.  As  the  Hamaker  constant 
increases,  the  particle  deposition  usually  is  enhanced.  It  should  be  noted  that  the  effects 
of  the  above  properties  is  much  more  significant  when  the  surface  interactions  play  more 
important  roles  in  the  particle  deposition,  or  when  the  height  of  the  total  interaction 
energy  barrier  is  higher  than  5 kBT.  If  the  deposition  is  convective  diffusion  limited, 
these  properties  will  not  have  any  obvious  effect  at  all.  Particle  deposition  becomes 
virtually  impossible  when  the  height  of  the  repulsive  energy  barrier  increases  beyond  20 
kBT.  The  effect  of  the  hydrodynamic  interaction  can  be  considered  as  a kind  of  repulsive 
“surface”  interaction.  It  exists  whenever  particles  move  close  to  the  collector  surface. 

Decreasing  particle  size  increases  particle  diffusivity  and  lowers  the  surface 
interaction  energy  barrier,  therefore  favors  deposition  of  colloidal  particles. 
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The  effect  of  the  operating  superficial  velocity  is  more  significant  when  the 
particle  deposition  is  convective  diffusion  limited.  The  particle  deposition  rate  increases 
as  the  superficial  velocity  is  increased. 


CHAPTER  7 

CONCLUSIONS  AND  SUGGESTIONS  FOR  FUTURE  WORK 

7-1  Conclusions 

In  this  research,  the  Finite  Element  Method  (FEM)  and  Effective  Medium 
Approximation  (EMA)  have  been  used  to  determine  the  fluid  flow  permeability  and 
particle  capture  efficiency  of  regular  array  and  random  arrays  of  cylindrical  and  spherical 
collectors.  The  FEM  takes  into  account  of  the  detailed  layout  of  the  packing  but  it  can 
not  handle  randomly  packed  collectors  yet  at  the  current  stage,  while  the  EMA  gives  an 
averaged  estimation  of  all  the  flow  parameters  but  is  not  able  to  provide  a resolution  to 
distinguish  the  detailed  structure  differences  and  orientations.  The  EMA  assumes  a 
model  system  in  which  a packing  element  (a  single  fiber  in  the  fibrous  medium  and  a 
single  sphere  in  the  granular  medium)  is  surrounded  by  a fluid  envelope  and  an  effective- 
medium  beyond  the  envelope.  It  integrates  the  important  features  of  both  the  cell  models 
and  Brinkman’s  model.  The  Stokes  equation  and  Brinkman’s  equation  are  solved  for  the 
fluid  envelope  and  effective  medium  regions,  respectively,  to  obtain  the  permeability  and 
close-to-surface  velocity  field  around  the  collectors.  The  convective  diffusion  equation  is 
then  solved  to  determine  the  particle  deposition  rate.  The  analytical  expressions  for  the 
permeability  and  particle  deposition  rate  are  derived  for  all  possible  cases  of  random 
packing  of  uniform  and  non-uniform  cylinders  and  spheres. 
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The  following  major  novel  conclusions  may  be  drawn  from  the  current  work: 

1 . The  numerical  simulation  using  the  finite  element  method  shows  that  not  only 
the  porosity  but  also  the  detailed  packing  structure  relative  to  the  flow  has 
significant  effect  on  particle  deposition  rate  in  a fibrous  medium.  The  more 
torturous  the  flow  streamline  is,  the  higher  the  particle  deposition  rate. 

2.  The  following  flow  permeability  equations  were  derived  based  on  the 
effective  medium  approximation: 


k = ~T 
a 


(1)  For  packing  of  uniform  circular  cylinders  parallel  to  mean  flow, 

a//RK,(a//R) 


-I1+— 
2 v 1- 


' 1 1 R2  R2 . A 

— r + InR 

4 4 2 


Voc// 


K0(a//R)+a//RK1(a//R)lnR 


(3-16) 


(2)  For  packing  of  uniform  circular  cylinders  perpendicular  to  mean  flow, 


ai=4(\-£) 


olrJ[-2oiR2  +3R4af  +16R2  -^]k,(^R)| 
Q {+8R'la1K0(a1R) 


(3-38) 


where 


Q = 


-a{  + 2a±R2  + a^R4  InR-aiR4  -ai  lnR  + 16R2a1  InR] 


RK,(a±R)+  4a2  R4  lnR  + 4aiRz  -3aiR4  -ai  +16R-  K0(a±R) 
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(3-35) 

(3)  For  plenary  packing  of  uniform  circular  cylinders  with  mean  flow 


parallel  to  packing  planes, 
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a2  =^(ai  +a//) 


(4)  For  completely  randomly  packed  uniform  circular  cylinders. 


2^2  1 2 

(X  — — Oli  H — 01// 

3 3 


(5)  For  randomly  packed  non-uniform  circular  cylinders. 


(3-48) 


(3-57) 
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If  the  size  distribution  is  not  very  broad,  the  following  equivalent 
radius  is  defined: 


ae  = 


Jo°  /o°  Za2Pa(a,Z)dadZ 


1/2 


v J0WJJ°  ZPa(a,Z)dadZ 


(3-71) 


All  the  equations  for  the  cases  of  uniform  cylinders  are  valid  for 
polydisperse  cylinders  if  we  replace  the  radius  of  the  cylinders  with  the 
ae  defined  above.. 

(6)  For  random  packing  of  uniform  spheres, 


a2  ..  18(1  - e) 
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and 


Q = 180R3  +24R5a2  -45R4a2  -9a2  +4a3  + 180aR4  +4R6a3 
+ 1 0(aR)3  -9Ra3  -9a3R5  +30(aR)2  - 1 80R3a 


(4-19) 
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(7)  For  random  packing  of  non-uniform  spheres, 

J0”aD(R.a)P(a)da  (4-26) 

k J0  a P(a)da 

The  thickness  of  the  clear  liquid  envelope,  8,  is  estimated  form  the 
overall  porosity  of  the  pack: 


oo 

Ja2 3Pa(a)da 

!_e  = _o (4-27) 

J(a  + 5)3Pa(a)da 

0 

R in  equation  (4-26)  is  obtained  by  R=(a+8)/a. 

3.  The  following  particle  deposition  rate  equations  were  derived  based  on  the 
Interaction  Force  Boundary  Layer  Approximation  (IFBLA)  for  surface 
interactions  and  Effective  Medium  Approximation  (EMA)  for  the  flow  fields: 
(1)  For  the  entrance  region  of  a random  packing  of  circular  cylinders 
parallel  to  mean  flow. 
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Ja 
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(5-9) 


(2)  For  random  packing  of  circular  cylinders  perpendicular  to  mean  flow, 


Sh  = 0.731A"3Pe1,3Tc(Pc) 


(1-38) 


where 
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Ac  = 
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(3)  For  plenary  and  completely  randomly  packed  circular  cylinders. 


Sh  = 0.731Aj./3Pe,/3yc(pc)x 


(5-15) 


where  x=0.8235  for  a plenary  random  packing  and  0.9107  for  a 
completely  random  packing. 

(4)  For  randomly  packed  non-uniform  circular  cylinders, 

J=O.731Xc0Ul':iZJ2,3/“A;,3[e,-^ja-2'3rc(e,a)Pa(a)PA(a)da 


(5-18) 

(5)  For  randomly  packed  uniform  spheres, 

Sh  = 0.653A'/3Pe1/3ys(Ps)  (1-38) 

where 
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Q -15a2R2  +12R6a3  +15a3R3  +3a2  +3a3R 

(6)  For  randomly  packed  non-uniform  spheres, 
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4.  In  a random  packing  of  circular  cylinders,  the  fluid  flow  relative  to  each  of  the 
cylinders  can  be  considered  to  be  the  vector  summation  of  a perpendicular 
component  and  a parallel  component  relative  to  the  cylinder  axis.  It  has  been 
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shown  that  the  contribution  of  the  parallel  component  to  the  overall  particle 
deposition  can  be  neglected  compared  with  the  perpendicular  component. 

5.  It  has  been  shown  that  if  the  height  of  the  repulsive  energy  barrier  in  the 
overall  surface  interaction  is  lower  than  5 kBT,  the  particle  deposition  can  be 
considered  to  be  convective  diffusion  limited  and  once  its  height  increases 
above  20  kBT,  the  particle  deposition  becomes  virtually  impossible. 


7-2  Suggestions  for  Future  Work 

7-2-1  Fluid  flow  through  fibrous  and  granular  media 
Packing  of  circular  cylinders  of  small  aspect  ratio 

In  a resent  work,  Rahli  et  al.  (1996)  found  that  the  flow  permeability  through  a 
random  packing  of  uniform  circular  cylinders  of  small  aspect  ratios  (e.g.,  short  fibers, 
length  / diameter  < 20)  can  be  20%  lower  than  that  of  a packing  of  long  circular 
cylinders.  It  is  believed  that  this  is  due  to  the  effect  of  the  fiber  ends,  which  have  been 
neglected  in  the  case  of  long  fibers.  If  the  total  surface  area  of  the  ends  is  comparable 
with  the  total  surface  area  of  the  fiber  sides,  theses  ends  will  have  a significant 
contribution  to  the  overall  particle  deposition  rate.  Theoretical  analysis  of  the  end  effect 
on  flow  permeability  and  mass  transfer  could  be  very  interesting. 
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Packing  of  non-spherical  collectors 

Granules  of  various  shapes  have  been  widely  used  in  industry  as  packing 
materials  for  catalysis,  absorption,  adsorption  as  well  as  filtration  processes.  There  have 
been  large  amount  of  flow  permeability  and  mass  transfer  data  available  for  packing  of 
these  granules  (Standish  and  McGregor,  1978).  The  theoretical  analysis  for  these  systems 
has  not  been  found  in  literature. 


7-2-2  Non-DLVO  interactions  and  colloidal  interactions  for  non-ideal  particles/collectors 
Fundamental  equations  for  non-DL  VO  interactions 

Until  now,  only  the  fundamental  theories  for  van  der  Waals  and  double  layer 
interactions  are  well  established.  Researches  on  the  other  surface  forces  such  as 
hydration,  polymer  bridging,  steric  interaction  and  hydrophobic  interaction  have  been 
very  active  in  the  recent  years.  The  ultimate  goals  of  these  researches  are  to  obtain  the 
fundamental  equations  for  calculating  these  forces. 

Colloidal  interactions  for  non-spherical  particles 

The  deposition  of  biological  particles  on  solid  surfaces  has  been  a very  hot 
research  topic  for  many  years  due  to  its  importance  in  health  related  applications 
(Dickinson,  1997).  The  interactions  between  these  non-regular  biological  particles  and 
solid  surfaces  ultimately  determines  whether  they  will  be  collected  on  the  surfaces. 
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Surface  interactions  for  rough  surfaces  and  surfaces  with  non-uniform  charges 

It  has  been  mentioned  in  Section  1-2-5  that  the  surface  roughness  and  non- 
uniformity of  charges  on  the  particle  and  / or  collector  surfaces  may  have  significant 
effect  on  the  particle  deposition  rate.  Although  there  has  been  some  work  done  on  this 
topic,  all  of  the  them  are  based  on  certain  assumed  patterns  of  roughness  and  charge 
distribution.  Research  is  needed  on  experimental  characterization  of  the  surface 
roughness  and  charge  distribution  of  the  actual  surfaces. 

Particle-collector  interactions  for  non-fresh  collectors 

As  particle  deposition  proceeds,  the  surface  characteristics  of  the  collectors  such 
as  the  surface  roughness  and  surface  potential  changes.  The  effect  of  these  deposited 
particles  on  the  surface  characteristics  and  the  fluid  flow  is  the  key  to  dynamic  modeling 
of  the  particle  deposition  process. 


APPENDIX 

INTRODUCTION  TO  SOLID  ANGLE 

A line  is  generated  by  sweeping  a point  linearly  through  space,  an  area  is 
generated  by  sweeping  a line  through  space,  and  an  angle  is  generated  by  sweeping  a line 
through  an  arc.  A solid  angle  is  generated  by  sweeping  two  “orthogonal1”  angles  (S  and 
<j) ) through  arcs,  as  shown  in  Figure  A. 

A solid  angle  is  two-dimensional  by  virtue  of  having  no  radial  position,  i.e.,  the 
solid  angle  remains  fixed  independent  of  radius.  This  provides  the  key  to  the 
mathematical  formulation  for  a solid  angle.  A differential  angle  is  defined  as  follows: 


d3  = — 


(Al) 


where  s is  the  differential  arc  length  and  r is  the  radius.  A differential  solid  angle  doo  is 
defined  by  scaling  up  the  differential  arc  length  to  a differential  area  and  dividing  by  r to 
make  the  solid  angle  independent  of  radial  position  (distance  from  the  origin)  as  follows 


dg)  = dA  = |BD|BC|  _ (rd9)(rsin  »#)  _ s|n  ^ 


(A2) 


Referring  to  Figure  A,  note  that  length  |BC|  is  generated  by  sweeping  line  |AB|=  r sinS 
(not  r),  through  angle  d<J>,  i.e.  the  revolution  is  of  a line  normal  to  the  vertical  axis  |AB|, 
about  the  vertical  axis  |OA|.  Contrariwise,  line  |BD|  is  generated  by  sweeping  |OB|=  r, 
through  angle  d9  about  the  origin  O. 
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Figure  A Solid  angle 

In  Figure  A,  9 is  the  “cone”  angle  and  <j)  is  the  “azimuthal”  angle.  While  angles 
are  measured  in  radians,  solid  angles  are  measured  in  steradians.  A hemisphere  contains 
2tc  steradians,  as  follows: 

A 0)  = f T sin  9d9dd>  = 2n 

J9= o Jb=o 


(A3) 
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